zasm 


1 49065 


Reports of the Department of Geodetic Science 
Report No. 198 


NASA OR- 
/ ^ / 7 ^/ 



Ivan h Mueller 
editor 


Prepared for 


National Aeronautics and Space Administration 
Johnson Space Center 
Houston, Texas 


Contract No. NAS 9-9695 
OSURF Project No. 2841 
Find Report 



f{N.ASA-CR-141751). , INVESTIGATIONS RELATED TO 
GEODEIXO CONTROL .ON THE MOON Final Report 
(Ohio -State Dni-v. .Research .Foundation) 

10.2 p 


The Ohio State University 
Research Foundation 
Columbus, Ohio 432 1 2 


N75-72999 


Unclas 

0.0/98 . • 09505 


V... * ■ 

< Cut 7/ * * 

REPRODUDED BY ~~ 

NATIONAL TECHNICAL 
INFORMATION SERVICE 

U. S. DEPARTMENT OF COMMERCE 
. . SPRINGFIELD, VA. 22161 




Reports of the Department of Geodetic Science 
Report No. 198 


INVESTIGATIONS RELATED TO GEODETIC CONTROL 
ON THE MOON 


Ivan I. Mueller 
editor 


Prepared for 

National Aeronautics and Space Administration 
Johnson Space Center 
Houston, Texas 


Contract No. NAS 9-9695 
OSURF Project No. 2841 
Final Report 


The Ohio State University 
Research Foundation 
Columbus, Ohio 43212 


February, 1973 



PREFACE 


This project is under the supervision of Professor Ivan I. Mueller, 
Department of Geodetic Science at The Ohio State University, Columbus, 
and is under the technical direction of Mr. Richard L. Nance, Code TF 541 
Mapping Sciences Branch, Earth Observation Division, NASA/OSC, Houston, 
Texas. The contract is administere d by the Facility and Laboratory 
Support Branch, Code BB 631 /B4, NASA/JSC, Houston, Texas. 


Preceding page blank in 



TABLE OF CONTENTS 


PREFACE iii 

TABLE OF CONTENTS v 

LIST OF TABLES vii 

LIST OF FIGURES viii 

1. INTRODUCTION 1 

2. DYNAMICS OF THE EARTH-MOON SYSTEM 4 

2.1 Orientation of the Earth by Numerical Integration 4 

F.A. Fajemirokun and Frank Hotter 

2.11 Introduction 4 

2.12 The Differential Equations of Motion of the Earth. . . 6 

2.13 Numerical Integration of the Differential 

Equations of Motion 11 

2.14 Numerical Experiments and Results 14 

2.141 Fitting the Numerically Integrated Eulerian 

Angles to Simulated Angles 15 

2.142 Comparison of Numerically Integrated Eulerian 

Angles with Classical Reductions 21 

2.2 Orientation of the Moon by Numerical Integration 27 

H.B. Papo 

2.21 Introduction 27 

2.22 General Principles 28 

2.23 Equations of Rotational Motion of the Moon 29 

2.24 Least Squares Analysis 32 

2.25 Comparison to an Analytical Solution 34 

2.3 Simulated Earth-Moon Environment 37 

H.B. Papo 

2.31 Introduction 37 

2.32 Qualitative Description of the Simulation 38 

2.321 The Earth 38 

2.322 The Moon 39 

2.323 Lunar Satellites 41 

2.324 Illuminating Source "Sun" 41 

2.325 General Notes 41 

Preceding page blank 


V 



2.33 Constants Used in the Simulation 42 

3. POSITION DETERMINATION FROM EARTH-MOON OBSERVATIONS 45 


3.1 lunar Far-Side Positions from Apollo Trans-Earth 

Trajectory Photography 

M. D. Sprague and W. Riotte 


3.11 Summary of Theoretical Approach. ...... 

3.12 Results of Simulated Block Triangulation Adjustments . 

3.13 Results from Apollo 15 Trans-Earth Photoqraphy . . . . 

3.14 Conclusions and Recommendations 

3.2 Lunar Near-Side Positions from Earth-Based Photography . . . 
H.B. Papo 

3.21 Introduction 

3.22 Basic Concepts 

3.23 Adjustment Theory 

3.24 Numerical Experiments 

3.25 Summary 


45 

47 

54 

59 

61 


61 

63 

67 

71 

73 


3.3 Positioning from Lunar Laser Ranging 
F. Fajemirokun and F. Hotter 


3.31 Introduction .... 

3.32 Numerical Experiments 

3.33 Results 


4. CURRENT WORK 


5. BIBLIOGRAPHY 


95 


5.1 Interim Reports Prepared Under the Contract 

5.2 Other References 


vi 



LIST OF TABLES 

Page 

2.1- 1 Gravitational Constants 13 

2.1- 2 Starting Initial Conditions for Three Test Cases 17 

2.1- 3 Recovered Initial Conditions for Three Test Cases 

After Two Iterations 18 

2.1- 4 Correlation Matrix for Adjusted Initial Values in 

the Simulated Case 18 

2.1- 5 Correlation Matrix for Numerically Integrated Eulerian 

Angle Fit to Classical Theory 24 

3.1- 1 Comparison of Standard Deviation for Five Representative 

NEC Points from Block Adjustments of the Simulated 
Apollo 14 Mission ' . 51 

3.1- 2 Standard Deviation of Selenographic Coordinates (in meters) 

for NEC Points Resulting from Sixteen-Photo Adjust- 
ments of Combined Mission Data 52 

3.1- 3 Summary of Results 12 Photo Block Adjustment (6 Iterations) 58 

3.2- 1 Covariance and Correlation Matrices for Solution of 

Network . I from 30 Bundles. Triangulation Points. ... 75 

3.2- 2 Solution Vector for Test (a) in Experiment (ii) 76 

3.3- 1 Geocentric Coordinates of Laser Stations 86 

3.3- 2 Selenodetic Coordinates of Lunar Retroref lectors 86 

3.3- 3 Variances of Parameters for Case 1, (One Laser Station 

Observing, Number and Period of Observations Varied) . 87 

3.3- 4 Variances of Parameters for Case 2, (One Laser Station 

Observing, Period of Observations Varied) 88 

3.3- 5 Variances of Parameters for Case 3 (Three Laser Stations 

Observing, Number and Period of Observations 

Varied) 89 



LIST OF FIGURES 




2.1- 1 The fundamental coordinate systems and the Eulerian 

angles 10 

2.1- 2 Residuals in 0 and i|> before and after adjustment 

(simulated case) 19 

2.1- 3 Residuals in $ before and after adjustment (simulated 

case) 20 

2.1- 4 Correlation between the Eulerian angles 23 

2.1- 5 Residuals of Eulerian angles • 25 


2.2- 1 Coordinate systems 

2.2- 2 Comparison to Eckhardt (1970) . . . . 

2.2- 3 Comparison to Eckhardt augmented for additive and 

planetary terms 


3.1- 1 Errors in X,Y,Z selenographic coordinates 

3.1- 2 Frame 2780 at TEI + 30 min 


3.2- 1 Earth-moon environment 

3.2- 2 

3.2- 3 

3.2- 4 Fundamental control network 

3.2- 5 From the earth the moon is seen through a narrow window . 


65 

68 

68 

71 

73 


viii 



1 . INTRODUCTION 

The main objective of this project was to design a system capable of 
providing geodetic information on the moon (and on the earth) consistent 
with the dynamics of the earth-moon system. This consistency is imperative 
since most modern observational systems between the earth and the moon 
(laser, differential VLBI, earth-based or orbital lunar photography, etc.) 
cannot be fully analyzed and utilized without taking into consideration all 
the factors which make any of them dependent on time, i.e., on dynamical 
factors. The most important group of factors on which the magnitude of the 
observations and their variations depend are the following: 

Group 1 . Parameters defining the position of the observation stations 
on the earth (or their motion) with respect to an earth-centered and earth- 
fixed coordinate system. 

Group 2 . Parameters defining the motion of the earth-centered and 
earth-fixed coordinate system with respect to an inertial frame of 
reference. • 

Group 3 . Parameters defining the motion of the moon-center with 
respect to an earth-centered and inertial ly oriented frame of reference. 

Group 4 . Parameters defining the motion of a moon-centered and moon 
fixed coordinate system with respect to an inertial frame of reference. 

Group 5 . Parameters defining the lunar stations (targets) or their 
motions with respect to the moon-centered and moon-fixed coordinate system. 
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In the classical geodetic/astronomic practice, the approximate number of 
parameters in the above groups would be the following: 

1. Earth-station parameters: 3 coordinates/station 

2. Parameters defining the orientation of the earth: 


— rotational velocity 1 

— precession (Newcomb) 9 

—nutation (Woolard) 109 

—polar motion ^ 50/year 


3. Parameters in an analytical lunar ephemeris: ^1 600 

4. Parameters defining the orientation of the moon (physical libration): 

^50 

5. Lunar station parameters: 3 coordinates/station 

It was clear from the beginning that it would be impractical to design 
a system which would attempt to correct the above several thousands of 
parameters simultaneously and at the same time possibly also try to 
recover certain systematic errors burdening the observations themselves. 

In order to reduce the number of unknown parameters and thereby make the 
approach more practical, it was decided in the beginning that the classical 
models for the orientations of the earth, the moon, and for the motion of 
the moon around the earth (groups 2-4 above) will each be replaced by 
sets of Eulerian angles and their time rates used in numerical integration 
procedures. Numerically integrated lunar ephemerides were already 
available (e.g., JPL) at that time. 

Another objective of the project is to find possible improvement in 
the earth-moon system described above when the different types of 
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observational systems are used in various combinations with each other 
rather than individually, which seems to be the current practice. In the 
beginning it was thought that this question could be best answered in a 
simulated and simplified earth-moon environment where the parameters in 
groups 2-4 can be adequately computed and thus the system can be 
controlled at will for an extensive simulation study. 

The project is not completed at the time of the termination of the 
present contract. It is continued under Contract No. NAS9-13093 
as OSURF Project No. 3487-A1. This report attempts to summarize those 
individual research items which have been completed and which will be used 
in the attempt to reach the original goals of the project. Most items have 
already been reported in more detail in the interim reports listed in the 
Bibliography. The report is in two main sections: The first (section 2) 
summarizes results related to the dynamics of the earth-moon system 
(groups 2-4), such as the simulated earth-moon environment, the numerical 
integration of earth and moon orientation angles; the second (section 3) 
summarizes some preliminary results in connection with position determina- 
tions on the earth and on the moon (groups 1 and 5). Section 4 is a summary 
of planned future work. Credit to individual contributors is given in each 
section, but it must be understood that in a university environment when 
the researchers are mostly graduate students and their advisors, the 
information flow is so free and the exchange of ideas so continuous that 
the names at the section headings are more indicative of responsibilities 
and not necessarily of individual achievements, though the latter obviously 
should not be ruled out either. 
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2. DYNAMICS OF THE EARTH-MOON SYSTEM 


2.1 Orientation of the Earth by Numerical Integration 
F.A. Fajemirokun^ and Frank Hotter 


2.11 Introduction 

The earth's rotational motion on its axis and its orbital revolution 
around the barycenter are of fundamental importance in astronomy. Both 
motions are complex and irregular due to the gravitational attractions of 
celestial bodies in the solar system, the tidal deformations of the earth, 
and the noncoincidence of the earth's axis of rotation with a principal axis 
of inertia. 

The changes in the direction of the earth's instantaneous axis of rota- 
tion in an inertial space have been obtained classically by dividing the 
variations into two parts— namely, the secular part (precession) and the 
periodic part (nutation). The differential equations of motion are 
analytically solved by the well-known method of "variation of parameters" 
and, through successive approximations, and after fitting actual observa- 
tions to the linearized mathematical structure, series expressions are 
obtained for the calculation of precessional and nutational elements. There 
are also variations in the angular rate of rotation (which affect the 
measurement and determination of time) and variations in the position of 
the instantaneous rotational axis within the earth (polar motion). 

‘Presently at the University of Lagos, Lagos, Nigeria 
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Presently, values of precession and nutation for any epoch calculated 
from the series expressions can be expected to be consistent (accuracy- 
wise) with present-day observation methods and instrumentation. Neverthe- 
less, the introduction of more precise instrumentation for future observa- 
tions, such as the laser and the Very Long Baseline Interferometry (VLBI), 
has recently reopened the possibility of more accurate determination of the 
parameters associated with the earth's rotational motion (such as precession, 
nutation and polar motion). Using these new types of observations, observa- 
tion equations expressing the observables as functions of parameters (which 
include those associated with the earth's rotation) can be used in an 
adjustment process to obtain corrections to existing values of the rotational 
parameters of the earth. 

The large number of parameters involved necessitates a large number of 

observations over an extended period of time. The coefficients of the 

series expressions for nutational parameters alone number 109, while the 

expressions for the precessional elements contain 9 coefficients. This 

large number of coefficients results from the analytical methods used in 

* 

solving the original equations of motion. Theoretically, the three 
second-order differential equations can be solved if six constants of 
integration are stipulated. 

In this paper, a method for solving the differential equations of 
motion of the earth is presented. The method is based on the numerical 
integration. The numerical integration yields the orientation angles of a 
coordinate system rigidly fixed to the earth's body with respect to a fixed 
celestial coordinate system. The parameters associated with the earth.' s 
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orientation in space could thus be reduced to six initial conditions, 
namely, the orientation angles and their time rates at a standard epoch. 


2.12 The Differential Equations of Motion of the Earth 

The motion of any rigid body under the influence of external forces can 
be represented as the resultant of a translation, with the velocity of the 
center of mass, and a rotation about an axis through the center of mass. 

The rotational motion is determined by the moments of the external forces 
about the center of mass and are well treated in standard textbooks on 
celestial and analytical mechanics. 

The force-moment equation for a rigid body can be written, according 
to Edwards [1964] as 


where 


-V ■£ 

L = 2 r. x F. 
n i i 


force momentum 


H = z m.r. xv. = rotational momentum 
n 1 1 i 

and 

m ... m are the masses of the particles which make up the rigid body 
i n 

r. ... r^ are the respective position vectors of the masses of particles 

from the center of mass of the rigid body 

F. ... F are the external forces acting on the particles 
i n 

V. V are the velocities of the particles referred to the fixed origin. 

1 n 


For the rotating earth, let U,V,W axes represent a Cartesian coordinate 
system fixed to the earth’s body and centered at the geocenter {for example, 
the 'average 1 terrestrial coordinate system). Further, let the X,Y,Z axes 
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represent an inertial coordinate system centered at the geocenter (for 
example, the mean ecliptic system of 1950.0). These coordinate systems are 
shown in Fig. 2.1-1. If OP is the, instantaneous axis of rotation, 
coinciding with neither the W nor the Z axes, then the earth's rotational 
velocity can be- represented by a vector W along OP. If i, j, k are unit 
vectors along the U,V,W axes respectively, then 

co = oj i + oj j + to k 2.1-2 

u v w 


Since 


and 


then 


v. = a) x r. 
i i 


-y -y -y 

r * ‘ u? + v? + w^ => u i + v j + w k 
i i j k 


H * H U 1 + H v J + H* k 

Differentiating H with respect to time: 


2.1 - 3- 
2.1 - 3 
2.1 - 4 


dH 

dt 


-y 

i 


dH, 

dt 


-5* 

di 

H — + 
dt 


dH, 

dt 


+ H, 


dj 

dt 


+ k 


dH, 

dt 


w 


dk 
H w dt 


2.1 - 5 


From Fajemirokun [1971] the above equation can be written in matrix notation 
as 

[H] = [Mj ] [»] =* [Mj] [«] 2.1 - 6 

where 


H 

u 

H v 

II 

1— t 

E 

A 

-F 

-F 

B 

-E 

-D 

9 

(0 = 

U u 

“v 

H w 


-E 

-D 

C 



“w 


and 

A,B,C,D,E,F are the six constants of inertia as defined in Edwards [1954]. 
Thus 
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4r- = M t w + to x H 
dt I 



If V represents the total potential energy of a system of n bodies of 
finite dimensions acting on the earth, and , X 2 , ... X n are the vectors 

■ 4 " 

from the geocenter to the centers of the bodies, then the force moment L is 
given as [Fajemirokun, 1971]: 

L = dV«| x X-, + dV 2 x X 2 + ... + dV n x X n 2.1 - 8 

where 




G is the gravitational constant, 

M.j is the mass of body i, 

R,. is the distance of body i from the geocenter. 
Consequently, 


h*i * *1 + ^- m i*2 * h + ••• + ir Vn * - 9 

1 2 


Equating equations 2.1 - 7 and 2-1 - 9, 


Mtui + a) x Mtoj “ -3G _J — MtXi x X-j + ... + — P— MjX n x X n 
1 1 L R5 1 R n J 


2.1 - 10 


which is the unsimplified dynamic equations of motion of the earth 
(rotational), assuming the earth is a rigid body. Equation 2.1 - 10 reduces 
to the well-known Euler's equations of motion [Plummer, 1960] if it is 


further assumed that the products of inertia (D, E, F) are zero. 

Through equation 2.1 - 10, the earth's rotational motion is represented as 
a rotation about an axis which constantly passes through its center of mass. 
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but whose position within the earth and. its direction in space varies from 
instant to instant. A geometrical equation can also be obtained by 
defining two coordinate systems, one fixed in space and the other fixed to 
the earth’s body. Then the motion of the earth-fixed coordinate system 
with respect to the space-fixed system' could be expressed mathematically. 

The orientation of the earth in space at any instant can be defined as 
usual by the three Eulerian angles [Edwards, 1964; Woolard, 1953] e, ip and 
$. These angles are shown in Fig. 2.1-1. From the figure, the following 
geometrical relationships follow: 

&> u = - § cos $ - sin 0 sin $ 

co = - 9 sin $ - ip sin e cos § -2.1-11 

co w = i cos 0 + $ 

or, in vector notation 

2 = S-j E 2.1-12 

where 



- cos $ 

- sin 

0 

sin $ 

0 


0 


to u 

S 1 ■ 

- sin $ 

sin 

0 

cos $ 

0 

, E = 

¥ 

, to - 

to v 


0 

cos 

9 


1 


$ 


“w 


The differentiation of equation 2.1 - 12 with respect to time gives 


* . t * 

co = S-jE + S-jE - S-jE + S 9 E 


2 2 


2.1 - 13 


where 


-cos 8 sin $ sin $ -sin 0 cos $ 
-cos 0 cos $ . cos $ sin 0 sin $ 
-sin e 0 0 


E 2 ■ 


0 Ip 

• « 

0 $ 

• • 
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p 



a) — Earth’s Rotational Velocity 
XYZ — Mean Ecliptic (1950. 0) Coordinate System 
UVW — Average Terrestial Coordinate System 
0 — Eulerian Angles 

Fig. 2.1-1 The fundamental coordinate systems and the Eulerian 
angles. 

Equation 2.1 - 13 expresses the earth's rotational motion geometrically in 
terms of the Eulerian angles and their time rates. Substituting these 
geometrical equations into the dynamical equations 2.1 - 10 and retaining 
E on the left-hand side, we obtain: 
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The above equation is a system of three second-order differential equations 
of motion, and the solution -for £ (the Eulerian angles) gives the orientation 
of the earth in space at any desired epoch. 


2.13 Numerical Integration of the Differential Equations of Motion 

In previous analytical solutions of the earth's differential equations 
of motion [Woolard, 1953], the earth's products of inertia (D, E, F) are 
assumed to be zero, and it is also assumed that the earth possesses rota- 
tional symmetry about the W-axis (i.e., A = B). Furthermore, the external 
forces acting on the earth are limited to those produced by the sun and the 
moon. - 

In this study, the moment of inertia (M^ ) matrix is taken to be a full 
one, and A is not necessarily equal to B. Also, in addition to the sun 
and moon, four other planets are used as sources of external forces in 
equation 2.1 - 14. These planets are Mercury, Venus, Mars and Jupiter. 

The three second-order differential equations of motion (equation 2.1 - 14) 
are transformed into. a form suitable for numerical integration by setting 
E 

Y = . 

E .• 

so that 

• • 

E E 

Y = 2.1-15 

E f(E,E) 

which represents a system of six first-order differential equations. These 
equations can be numerically integrated by any available numerical 
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integration algorithm, based on any standard numerical integration method. 
The method used in this study is the variable-order Adams predictor- 
corrector method, developed and used by the Jet Propulsion Laboratory (JPL) 
in their Orbital Determination Programs. A description of the integration 
algorithm and its special advantages are contained in [Krogh, 1969]. 

In order to obtain the earth's orientation angles at a given epoch t, 
the six first-order differential equations of motion (represented by 
equation 2.1 - 15) are numerically integrated from a starting epoch t Q . The 
six constants of integration are the starting values: 



which may also be referred to as the initial conditions. e 0 , tJ> 0 and $o are 
calculated (for epoch to) using equations derived in Fajemirokun [1971], and 
their time rates are obtained through numerical differentiation. Also, in 

-3- 

order to numerically integrate equation 2.1 - 15, E must be evaluated from 
equation 2.1-14. E is dependent on certain parameters in addition to the 
values of 0 , if), and $. These parameters include the earth's moment of 
inertia matrix, the gravitational constants of the six celestial bodies, and 
the geocentric position of these celestial bodies. 

Using the coefficients of spherical harmonics taken from Rapp [1969J, 
the moment of inertia matrix used in these studies is: 


0.32967 

-0.1692 * 10" 5 

0 

0.1692 x 10 

0.32968 

0 

0 

0 

0.33076 
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where 


1 

K ~ Ma* 
e 

a g , M are the earth's semi -major axis and mass respectively. The geocentric 
position of the sun, moon and. the four planets were taken from the JPL 
Development Ephemeris 69 (DE-69) tape [O'Handley et al., 1969] which is 
obtained through numerical integration and is believed to be gravitationally 
consistent. The planetary masses used are from the JPL system of planetary 
masses [Melbourne et al., 1968] and Table 2.1-1 gives the gravitational 
constant for each astronomical body. 


Table 2.1-1 

Gravitational Constants 


Body 

Gravitational Constants 
(Km 3 / day” 2 ) 

Earth 

0.2975542 x 10 ls 

Moon 

0. 3659906 x 10 14 

Sun 

0.9906936 x 10 31 

Mercury 

0. 1655848 x 10 lS 

Venus 

0. 2425068 x 10 ls 

Mars 

0. 3197127 x 10 16 

Jupiter 

0. 9458682 x 10 18 


The values obtained for the integrated angles and their time rates at - 
any desired epoch depends, to some extent, on values chosen for the initial 
conditions and the parameters of integration (e.g., inertia matrix, gravita- 
tional constants, etc.). It is therefore necessary to provide a means 
through which the estimated values of the initial condition and other 
parameters can be corrected in an adjustment process involving some form of 
observation. Such a procedure is described in Fajemirokun [1971] and utilizes 
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state transition and parameter sensitivity matrices. The matrix of partial 
derivatives 


3E 

3E 

3E q 

3E 0 

ii 

3E 

a? 

m 

o 

3Eo 


is the state transition matrix, while the parameter sensitivity matrix is 
given by 

3E 

3a 

V = . 

3E 

3a 

where a represents the list of parameters. 


2.14 Numerical Experiments and Results 

In the previous sections of this study, it has been proposed to calcu- 
late the earth's orientation in space by numerically integrating the rota- 
tional differential equations of motion of the earth. The equations of 
motion were derived with only one assumption — that of the earth's rigidity. 
Such a new practical approach to an old problem as proposed here demands 
some confirmation as to the validity of the method numerically, and some 
verification of the correctness of the derived equations. Consequently, 
some numerical experiments were performed mainly to verify the equations, 
the computer programs and the procedure for adjusting the initial conditions. 

An important aid in verifying the computer programs was the Simulated 
Earth-Moon Environment Data (hereafter referred to as the simulated data) 
which was created by Papo £1971]. The generated data consists of the 
geocentric position and velocity of the moon, the orientation angles of 
the earth and the moon and their time rates, the state transition matrices 
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of both the earth and the moon and their parameter sensitivity matrices. 
These quantities were obtained by numerical integration of equations based 
on a moderately complex model of the earth-moon dynamic system consisting 
of a rotationally symnetric rigid earth and a perfectly rigid moon whose 
dynamic shape is that of a tri axial ellipsoid. Other details of the mathe- 
matical formulation of the equations of motion of this simplified earth- 
moon system are contained in the work of Papo [1971]. 

A simulated ephemeris of the earth and the moon was created for a 
period of one year beginning at 2440222.5 JD (1969.0). The ephemeris 
contains numerically integrated geocentric positions and velocities of 
the moon and the Eulerian angles and their time rates for the earth and 
the moon, recorded at half-day intervals. In addition, a fifth-order 
modified Everett interpolation formula [O'Handley et al . , 1969] was avail- 
able for use in interpolating the quantities at epochs which fall between 
tabulated values. 

The following two groups of numerical experiments are reported in this 
section: 

(a) Fitting the numerically integrated earth's Eulerian angles to those 
obtained from the simulated data. 

(b) Comparing the numerically integrated Eulerian angles of the earth to 
their counterparts obtained through classical methods. 

2.141 Fitting the Numerically Integrated Eulerian Angles to 
Simulated Angles. 

This experiment was performed in order to provide an independent check 
on the equations of motion of the earth as developed in this study and the 
computer programs based on them. Since the simulated data was used, the 



computer program written for the real case was slightly modified to 
accommodate the following restrictions imposed by the simulated environment 
model : 

(1) In the moment of inertia matrix, A - B, and D = E = F = 0. 

(2) Only the moon is the external celestial body whose potential affects 
the rotation of the earth. 

(3) The geocentric position of the selenocenter at any epoch is that given 
by the simulated data rather than that obtainable from a real lunar 
ephemeris. 

The equations of motion were integrated using the DVDQ subroutine [Krogh, 
1969]. 

A comparison of integrated angles with their counterparts from the 
simulated data shows differences which were less than 0.0001 for 0 and 
O'.'OOl for iji and $ over the one year period. These results indicate perfect 
agreement between the two numerical integration programs (one for this study 
and the other used in creating the simulated data) as well as the equations 
on which the computer programs were based. 

The adjustment program and the ability of the adjustment method to 
recover the initial conditions were also tested numerically. For this pur- 
pose, the theoretical initial conditions (i.e., values of the integrated 
quantities at the initial epoch of 2440222.5 JD), read from the simulated 
data, were varied. Three test runs were made in which 

(1) only the initial angles were varied, 

(2) both the initial angles and their time rates were varied, 

(3) only the time rates of the angles were varied. 

The adjustment was performed over an interval of forty days beginning at 
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2440222.5 OD and ending at 2440262.5 JD. Values of the angles integrated 
with wrong initial conditions were compared with the "true" values at half 
daily intervals. Table 2.1-2 shows the "correct" initial conditions and the . 
initial conditions used in each of the three cases listed above. 


Table 2.1-2 

Starting Initial Conditions for Three Test Cases 


Symbol 

Correct 
Initial Values 

Initial Values 
for Case 1 

Initial Values 
for Case 2 

Initial Values 
for Case 3 

0 

0.4091596226 

0.4091693189 

0. 4091644707 

0.4091596226 

* 

-0. 131898 xlO 5 

-0. 110152 x Iff 4 

-0. 616711 xlO 5 

-0.131898 x10 s 

§ 

4. 8950936587 

4. 8951033498 

4. 8950985068 

4.8950936587 

0 

0. 659363 xlO 7 

0. 659363 xlO 7 

0.725299 xlO 7 

0. 725299 xlO 7 

i 

-0. 300932 x10 s 

-0.30093 xlO 3 

-0. 303942 x 10 s 

-0.303942 x10 s 

4 

6. 3003883741 

6. 3003883741 

6. 3003943826 

6. 3003943826 


(0, i/>, $ in radians, 9, 4> in rad/day) 


Table 2.1-3 presents the adjusted values of the initial conditions 
for the three cases after two iterations. Comparison of the recovered 
parameters with their "correct" values shows agreement up to the ninth 
decimal places (equivalent to 0'.'0002 for the angles and 0!!0002/day for 
their time rates). The residuals obtained after the adjustment (by 
comparing integrated angles with counterparts from the simulated data) were, 
for all three angles and in all three cases, less than O'.'OOOl throughout 
the adjustment interval. Figs. 2.1-2 and 2.1-3 show the residuals after 
adjustment for the most general case (case 2). 
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Table 2.1-3 


Recovered Initial Conditions for Three Test 
Cases After Two Iterations 


Symbol 

Correct 
Initial Values 

Recovered 
Initial Values 
for Case 1 

Recovered 
Initial Values 
for Case 2 

Recovered 
Initial Values 
for Case 3 

e 

0.4091596226 

0.4091596227 

0.4091596227 

0. 4091596226 


-0. 131898 x10 s 

-0. 131869 x 10 s 

-0. 131880 x10 s 

-0.131860 xlO 5 

<£> 

4. 8950936587 

4. 8950936584 

4. 8950936585 

4.8950936584 

• 

e 

0. 659363 xlO 7 

0. 652045 xlO 7 

0. 654780 xlO 7 

0.649609 xlO 7 

i 

-0. 300932 x10 s 

-0,298465 xlO s 

-0.299432 x10 s 

-0.300307 

& 

6. 3003883741 

6. 3003883718 

6. 3003883727 

6. 3003883735 


(Q, ip, <£ in radians, 9, ip, 4> in rad/day) 


As can be seen from Table 2.1-4 which presents the correlation matrix, 
a few of the parameters in the solution have poor separation from others. 

In particular, there is a strong correlation between $ and 4 and moderately 
high correlation between 9 and if). In spite of these high correlations, the 
result of this experiment shows that in the adjustment performed, the 


Table 2.1-4 

Correlation Matrix for Adjusted Initial Values 
in the Simulated Case 
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Fig. 2.1-2 Residual's in 0 and i|» before and after adjustment 
( simulated case). 




Fig. 2.1-3 Residuals in $ before and after adjustment 
(simulated case). 




solution for all the parameters actually converges to their true values. 
This particular adjustment performed over forty days with only two itera- 
tions took approximately five minutes to run on an IBM 360/75 computer. 


2.142 Comparison of Numerically Integrated Eulerian Angles with 
Classical Reductions. 

In order to further test the numerical theory, experiments were 
conducted to compare the integrated Eulerian angles with similar angles 
computed from classical theory. Three problems became evident: 

(1) First, the integrated angles refer to the body axis of the 
fictitious solid earth. Hence, to compare with classical precession 
nutation theory, it is necessary to refer the angles to the rotation axis. 
If we define the rotations to the rotation axis by 


tan"^ 

r 

CiJ'j 


l“3J 

tan - ^ 

* ^ 


refer 

to 


2.1 - 16 


, o> , a) given by equation 2.1 - 11 The 

I ' L' 6 U V Z 

observed celestial coordinate system is then given in terms of the 
Eulerian angles by 


- R-j { — 3 ) R2 ( ot ) ^3 ( ) Rq (~o)R3 


OBS 


2.1 - 17 


1950.0 


Correspondingly, the same coordinates are given through classical reduction 
procedures using the Greenwich Apparent Sidereal Time (GAST), the precession 
and nutation matrices (P, N) and the obliquity of the ecliptic (e Q ) at the 
epoch 1950.0 as {Mueller, 1969] 
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In order to compare the two theories, it was assumed that the 3x3 matrix 
products from equations 2.1 - 17 and 2.1 - 18 were composed of three angles 
between the mean ecliptic system of 1950.0 and the observed system as 

= R 3 (GAST)N P R^-h) 

u » I “ 1 

^3 (4 s )Ri (— ©) R 3 = R] (~B ) R 2 C 01 )^3 )R"] (-©)R 3 (4».) 

Then the differences between the angles were minimized by least squares 
theory to adjust the desired initial state vector and/or parameters in the 
numerical theory. 

(2) Second, the gravitational constants and moment of inertia ratio 

(C - A)/C in the conventional theory are defined by the constant of nutation 
of 9'.'210 at 1900 and the constant of general precession in longitude 
5025V64 per tropical century at 1900 [Woolard, 1953, p. 124]. Therefore, 
to make the comparison realistic these values were used in the numerical 
integration and as in Woolard' s theory the principal moments A and B were 
assumed to be equal , with nondiagonal elements (D, E, F) assumed to be 
zero. 

(3) The third difficulty arose from assuming A = B and D = E = F = 0. 
From equation 2.1 - 10 it can be shown that 

= 0 


or 


CO- 


= constant = k 


2.1 - 20 


Hence, from equation 2.1 - 11 we have 

tog = k = ^ cos 0 + <f> 2.1-21 

Therefore, if we "set" an initial state vector defining ifi, 0 , and <f>, the 
constant is determined and we would expect ij>, ^ to be highly correlated as 
shown in Fig. 2.1-4. To overcome this problem <j> was assumed to be dependent 
on \ p and only, the values 0, 4>, ^ in the state vector were adjusted. 

Further, it was assumed that over the period of integration t 0 to t: 
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SECONDS OF ARC DESCREPANCIES 



INTERVAL (JULIAN DAYS) 


Fig. 2.1-4 Correlation between the Eulerian angles. 



2.1 - 22 


$ = 



so that 


k(t - t 0 ) 


94 ) 

9k 


t - to 


r^. 

ijt cos 0 dt 

. 

to 


2.1 - 23 


With the above factors introduced in the numerical theory, an integra- 
tion/adjustment was performed over a 100-day period with comparisons to 
classical theory as outlined in equation 2.1- 19 made at 0.25 day intervals, 
the sum squares of the residuals prior to adjustment was 0.3263 seconds 
squared and after adjustment 0.0001 seconds squared. The correlation 
matrix is given in Table 2.1-5 and a plot of the residuals in Fig. 2.1-5. 


Table 2.1-5 

Correlation Matrix for Numerically Integrated Eulerian 
Angle Fit to Classical Theory 
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SECONDS OF ARC DESCREPANCiES 



INTERVAL (JULIAN DAYS) 


Fig. «!.l-5 Residuals of Eulerian angles. 



2.2 Orientation of the Moon by Numerical Integration 


H.B. Papo 


1 


2.21 Introduction 

The establishment of fundamental control on the moon is considered 
inseparable from the lunar epheneris and the orientation of the moon in 
space. In the ideal case the control network, the lunar ephemeris and 
the orientation of the moon should form a consistent set of parameters and 
data. Analysis of the available lunar ephemerides led us to the adoption of 
JPL's DE-69 ephemeris as a basis of our system. For the orientation of the 
moon a new solution was developed as an integral part of a comprehensive 
data reduction system. As is well known, the model of physical librations 
of the moon together with Cassini's three laws provide the tool for defining 
the orientation of the moon in inertial space. 

The problems associated with conventional analytical solutions for 
the physical librations of the moon which we addressed in our study are 
as follows: 

(1) Difficulties introduced by the linearization of Euler's dynamical 
equations. 

(2) Insuring an absolute consistency between the physical libration model 
and any ephemeris of the moon and the sun— in particular, a numerical 
ephemeris like JPL's DE-69. 

(3) The need to identify and select a set of parameters which can be 
^Presently at the Technion, Israel. 

Preceding page blank 
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determined in a least squares analysis of given observational material. 
Our present solution has succeeded in providing a. satisfactory answer 
to the foregoing problems. However, we have made no provisions as yet to 
account for the partially elastic behavior of the moon. 

2.22 General Principles 

As a model for developing our approach, we considered procedures 

analogous to orbit determination by numerical integration. The simplest 

definition of such a model is as follows: 

The physical libration angles and their time rates (a total 
of six quantities) at an arbitrary epoch are assumed to be known. 

Using those as initial values, the physical libration angles at. 
another epoch are obtained by numerically integrating the equations 
of rotational motion of the moon. 

Thus, a solution for the physical librations of the moon by numerical 
integration comprises the following: 

(1) Derivation of differential equations of rotational motion of the 
moon. 

(2) Development of computational procedures for the solution of a set of 
initial values of the physical libration angles and their time rates 
at some epoch, as well as the solution of a number of physical 
constants which govern the equations of rotational motion of the moon. 

The coordinate systems involved in our analysis are the following 
(Fig. 2.2-1): The origins of both XYZ and xyz Cartesian coordinate systems 

is set at the center of mass of the moon. The XYZ system is parallel to 
the mean ecliptic system: the X axis points to the mean vernal equinox of 
date while the Z axis is normal to the plane of the ecliptic. The xyz 
system, to be referred to as the selenodetic system, is fixed with respect 
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z 



Fig. 2.2-1 Coordinate systems. 

to the body of the moon and its axes x,y,z coincide with the principal axes 
of inertia of the moon. The transformation between the two systems is 
carried out through the three Eulerian angles <j>, if», 0 . The physical 
libration angles t, a, p, are related to the three Eulerian angles <j>, ty, 0 
(equation 2.2 - 1) through the respective mean longitude of the moon (L), 
the longitude of ascending node of its mean orbit (ft) and through the 
mean inclination of the lunar equator with respect to the ecliptic (I): 


T 


$ + 


L - TT 

a 

- 


- 

ft 

P 


0 


I 


2.23 Equations of Rotational Motion of the Moon 

The rotational motion of the moon, regarded as a rigid body, around 
its center of mass and with respect to an inertial frame of reference is 
described by Euler's dynamical equations: 
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2.2 - 2 


"x 


a 0 0 
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' 

yz 


yz 


«f 

£ 

N 

U) y 

= 
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xy 
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X 

3 


The rotational velocity components of the moon around the x,y,z axes 
are denoted as a> x , Wy, t« 2 , while u x , a> y) w z are the respective time rates 
of those velocities. 

The coefficients a, B,y are ratios between the principal moments of 
inertia of the moon as follows: 


- C - B 
a A 


B = 


C - A 


_ ' B - A 


where A, B, C are the moments of inertia about x, y, z respectively. 

The gravitational constants of the earth and of the sun are given as 
2 2 

k E and k S respectively. 

The selenodetic coordinates of the centers of mass of the earth and 

of the sun are given as [x y z]_ and [x y z] while r and r are the 

E 5 E b 

respective distances moon-earth (r £ ) and moon-sun (r^). 

The relationship between the rotational velocities ai x , a z and 
the Eulerian angles and their time rates <j>, i|>> 0, 9 ifij e is given by 
Euler's geometric equations: 


“x 



* 

-0 


e x 

) 

0 

w y 

- Rg('i') R-|{“0) 


0 

+ Rg(<l») 

e y 

l + 

0 

h 



i 


e z 

) 



2.2 - 3 
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The symbols R 3 (<j>), Rj (- 0 ) 9 etc. represent orthogonal rotation matrices 
in the conventional sense. 

The only new feature in these equations is the vector [e , e v , e_] 

x y z 

which represents the small rotational velocity of the X,Y,Z system with 
respect to an inertial frame of reference. This motion results from the 
motion of the ecliptic and the regression of the vernal equinox. We 
should note that e x and e y are relatively small quantities. 

We differentiate Euler's geometric equations with respect to time 
and substitute the results in Euler’s dynamical equations. After some 
regrouping, the results are the second-order differential equations of 
the Eulerian angles of the moon: 


= W 


0 

•e 

0 


,3k 2 E 

yT5 


yz 

xz 

xy 


+ 3k 2 S 


yz 

xz 

xy 


“y “ z 


“x “z 


0) (J) 

x y 



-0 


e 

_ -J 



x / 

3R3 (<j>) Ri(-e)<fr + R 3 U) 3Ri (- 0 ) 0 ' 

0 

+ R 3 (ip) 

e / 

3<J> 36 

* 


y l 

1 o 1/ 



e z ) 


R 3 ($) Ri(“Q) 3Ra(4>) 


2.2 - 4 


All the symbols in these equations have been explained already except for 
W which is a 3 x 3 matrix composed of functions of 8 and <f>. 

In principle, these three equations can be used for numerical integra- 
tion of the Eulerian angles. However, due to the large numerical values 
of <f>, \p, e and also due to a strong negative correlation between 4) and \p. 


31 



it is advantageous to perform another transformation and obtain differential 
equations of the physical librations. This is similar to a transformation 
from Cowell's to Encke's equations of motion of a satellite. 

By repeated differentiation of equations 2.2 - 1, we arrive at the 
final result to be referred to as the differential equations of rotational 
moti on : 


T 


• • 

<{, + !{, 


L 

a 

= 

II 

* 

- 

ft 

P 


0 


0 


where L and ft are composed of power series of time with fixed coefficients. 
2.24 Least Squares Analysis 

The numerical solution of the differential equations of rotational 
motion of the moon (2.2 - 5) depends on two sets of parameters which are 
subject to estimation, namely: 

(1) The values of the physical libration angles and their time rates at 
an arbitrary epoch. Those values will be referred to as the initial 
values (? Q ). 

(2) A number of constants which appear on the right-hand side of the 
differential equations. 

We should note at this point that the ephemeris data and also the 
gravitational constants of the earth and the sun which are associated v/ith 
the ephemeris are treated in our analysis as being fixed. The same applies 
to the vector [e , e e l. 

A J 

The only parameters which are subject to change and which are being 
estimated in the least squares sense are the ratios a, 3, y. As the three 
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are interrelated by a well-known identity, we consider only two quantities 
to be referred to as the physical constants (x). Those could be a and 3 
or the traditional f and y where 

f - SL 
6 

The above parameters can be estimated from the analysis of any type of 
observational material which is, in general, sensitive to the orientation 
of the moon. Heliometer observations, photography by long focus astronomical 
telescopes or laser ranging to retroreflectors placed on the moon surface 
could be used for this purpose. 

Each type of observation is associated with a typical set of additional 
parameters which pertain to the instrumentation and to the particular 
observational mode being used. This last group of parameters will be 
referred to as observation constants (y). 

Each observation can be modeled by the generally nonlinear model F: 


F = F(c!, * a , u a ) 


2.2 


According to the standard procedure in least squares estimation, we proceed 
by linearization of the model F in terms of the parameters and obtain the 
observation equations: 


w 3F „ _l 3f , j. 3h . t- fr-0 ,0 o N r 

V = ^T?o + -57a * + IwT v + F(£ , * » v ) - F 


dF 


8F 


3 ^ 


2.2 - 7 


3X a 3y a 

We will pay closer attention to the coefficients of S 0 and X. The 
partial derivative 3F/3y does not present any particular problem, and it 
is irrelevant to the numerical integration process. The two coefficients 


9F/3? 0 and 3F/8X can be broken into the following forms: 
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3F 

K 


2.2 - 8 


. - lE . = A U 

9? 9?J £ 

3F _ 3F . 3£ _ a p 

K* ~ 3? " V 

3F/3? are the partial derivatives of the observational model F with respect 
to £ — the current values of the physical libration angles and their time 
rates. These partial derivatives (3F/3£), like the ones above (3F/3y), are 
irrelevant to the numerical integration procedure. 

The two matrices U and Q represent the partial derivatives of the 
current physical libration angles and time rates (?) with respect to the 
initial values fe ) and again with respect to the physical constants (x), 
respectively. In the literature on orbit determination these two matrices 
are referred to as the state transition (U) and the parameter sensitivity 
(Q) matrices. We should notice the important fact that U and Q are inde- 
pendent of the particular type of observation being used. 


2.25 Comparison to an Analytical Solution 

To test our solution we chose as a standard Eckhardt's 1970 model for 
the physical librations as published in [Eckhardt, 1970, p. 267]. After an 
extensive analysis we selected the function to be minimized by the least 
squares procedure as the sum of squares of Ax, a a sins and Ap. Here the 
a's represent the differences between the numerically integrated and 
Eckhardt's values at each epoch of comparison. 

We ran an adjustment over about 1400 days, taking Eckhardt's angles 
at intervals of three days. The residuals after the adjustment, reduced 
to an amplitude of 2 to 3 seconds of arc, except for the residuals in t 
which displayed a clearly periodic character with large amplitude and a 
period of roughly three years. The adjusted values for 6 and y agreed 



with Eckhardt's values as given in the paper mentioned above. 

We changed the initial epoch by several hundred days, adjusted again, 
and found that the residuals persisted with only minor differences. This 
demonstrated the stability of our solution (see Fig. 2.2-2). 



It was thought that the large differences in t were due to some inconsistency. 
The integrated physical libration angles were obtained through the use of 
the JPL DE-69 ephemeris, while Eckhardt used Brown's ephemeris without the 
additive and the planetary terms in the lunar theory [Eckhardt, 1967, p. 51]. 
Eckhardt's angles augmented by the inclusion of the effect of the additive 
and planetary terms [Williams, 1973] were used again as a new standard. 

The new adjustment over the same interval as above resulted in a 
significantly improved fit. The residuals in At and in Ap dropped below 
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Fig. 2.2-3 Comparison to Eckhardt augmented for additive and 
planetary terms. 

We concluded our experiments of fitting to an analytical solution, 
having shown that the two methods are fairly compatible over an extended 
period of time. The next test of our method would be made through the 
analysis of real data which is sensitive to the physical librations of 
the moon. 
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2.3 Simulated Earth-Moon Environment 


H.B. Papo^ 


2.31 Introduction 

As in any simulation, the primary goal is to generate absolute data 
which will be entirely self-consistent and, in addition, will resemble as 
closely as possible the true environment. 

After exploring briefly several possibilities the way chosen was to 
use numerical integration of the equations of motion (translation as well 
as rotation) of the earth, the moon, and a satellite orbiting the moon. 

The following sets of data are generated: 

(a) Geocentric ephemeris of the moon in a fixed celestial coordinate 
system. 

(b) Eulerian angles for orienting an earth-fixed coordinate system vs. 
the same inertial celestial system mentioned in (a). 

(c) Eulerian angles for orienting a moon-fixed coordinate system to the 
inertial celestial system. 

(d) Selenocentric ephemeris of a moon-satellite in the inertial celestial 
system. 

(e) Range and range- rate from four tracking stations on earth to the 
satellite of the moon. 

■j 

Presently at the Technion, Haifa, Israel. 
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(f) Optical bundles of rays from three observatories on earth or from a 
satellite orbiting the moon. Each ray represents the unit vector of 
the line connecting the projection center to a given triangulation 
point on the moon in a specially oriented reference frame. The case 
is demonstrated in the diagram below. 

PC - projection center 

B 1 ,B 2 ,E 5 - reference frame oriented vs. 

the inertial coordinate system 
through three Eulerian angles 

TP - tri angulation point on the 

moon, M 

axis - points to the moon's center 

B 2 axis - is parallel to the XY plane of 
the inertial system 

B 3 axis - completes a right-handed 
Cartesian system 

2.32 Qualitative Description of the Simulation 

The only celestial bodies considered in the system are the earth and 
the moon. The fixed (inertial) celestial coordinate system is chosen 
parallel to the mean ecliptic system of a given standard epoch. However, 
no stars are provided to serve as reference for inertial space. 

2.321 The Earth. 

The earth is a perfectly rigid body having a topography identical to 
the true earth. A Cartesian coordinate system centered at the earth mass 
center and coinciding with the principal axes of inertia is defined fixed 
to the earth body and following identically the UVW — average terrestrial 
system of the true earth: W passing through the CIO pole, plane UW 
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containing Greenwich. The equatorial semi -axis and the flattening of the 
reference ellipsoid are the same as the ones adopted by IAU in 1964. 

Four radio tracking stations (range and range-rate) and three optical 
observatories are defined at locations where actual and presently active 
stations are engaged in observations of the' appropriate type. 

Dynamically the earth is rotationally symmetric. The spherical 
harmonics expansion of the earth gravity field has only one coefficient- 
02 .--. different from zero. The mass of the earth and the value of J 2 are 
identical with the values adopted by the IAU in 1964. Choosing a value for 
the ratio between the polar and the equatorial moment of inertia (C - A)/C 
was a rather complicated matter. The value adopted finally was the one 
consistent with the nutation constant in obliquity N = 9V21. The value 
calculated by using Sir H. Jeffrey's formulae differs slightly from the 
one adopted. 

The orientation of the UVW system was taken from the American Ephemeris 
and Nautical Almanac series for the sun. Its rate of diurnal rotation was 
determined in accordance with Newcomb's equation. These served as the 
initial values in the integration. 

2.322 The Moon. 

The moon is a perfectly rigid body. Its topography is that of a sphere 
with radius equal to the mean radius of the true moon. Its center of mass 
is shifted with respect to the center of the sphere towards the earth. On 
the "front" side of the moon (facing the earth) 30 triangulation points are 
chosen evenly spaced and in areas on the true moon where there is an 
abundance of crater! ets of the 3-7 km diameter range. The crater Mosting A 


39 



and also the crater Bruce in Sinus Medii are among the 30 points. A 
Cartesian coordinate system fixed to the body of the moon, centered at its 
mass center and following the orientation of the selenographic coordinate 
system (x,y,z) is defined. The xyz axes are identical With the principal 
axes of inertia. 

The mass of the moon is related to that of the earth by the ratio 
1/m identical to the one adopted by IAU in 1964. 

Dynamically the moon is a more complicated body than the earth. It is 
a tri axial ellipsoid with a set of 12 point masses superimposed on its 
triaxial dynamic shape. In terms of coefficients of spherical harmonics 
expansion C 2 ,o and C 2 ,2 ane the only nonzero second-order coefficients. 
Their values were taken identical to the ones used by NASA. The mascons 
are selected in general at locations corresponding to actual mascons as 
published recently. There are several additional mascons introduced, some 
of them negative, selected to satisfy certain conditions. The conditions 
that the total of 12 mascons satisfy are that their total mass and also 
their first and second moments are all zero. This was found necessary in 
order to retain a basically rough gravity field as far as a close moon 
satellite is concerned and at the same time not to complicate unnecessarily 
the equations of motion of the moon vs. the earth. 

The three moments of inertia (A,B,C) are consistent with the C 2j0 
C 2 , 2 values and correspond to a ratio e = (C - A)/8 identical to the 
presently accepted value of 6. 

The initial orientation and rates of rotation of the moon's xyz system 
are taken from Eckhardt's model. 
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2.323 Lunar Satellites. 

A spacecraft of negligible mass is defined in which no parallax exists 
between its mass center, the calibrated point of the transponder for radio 
measurements and the principal point of its camera. No orientation jets or 
other physical effects disturb its perfectly gravitational motion. In 
designing the various satellite orbits no consideration is given to their 
trajectories prior to the circumlunar orbit. 

2.324 Illuminating Source "Sun". 

An illuminating nongravitational "sun" is defined rotating on the XY 
fixed celestial plane at an infinite distance from the earth-moon system. 

Its angular rotation speed corresponds to that of the mean sun as given by 
Newcomb's equation. 

2.325 General Notes. 

(a) In generating the optical, range and range-rate data all the 
practical considerations of visibility, illumination of target, photography 
away from the sun, minimum altitude angles for earth photos, maximum field 
angle of camera, etc. are enforced. If any of the conditions outlined 
necessary for the existence of an observation is violated, the tentative 
observation is replaced by zero. 

(b) The computer programs developed to generate the simulated environ- 
ment are based on Encke-type integration, i.e., only the perturbations to 
the reference cases are integrated, thus increasing the number of signifi- 
cant figures in the quantities needed in the simulation. 
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2.33 Constants Used in the Simulation 


Earth 

K 2 E = .297556 x TO 16 km 3 /day 2 

J 2 = .0010827 

a = 6378.16 km 

f = 1/298.25 

.00327802 

to = 6.300388098 rad/day 


gravitational constant 

second zonal harmonic 

equatorial .semi -diameter 

flattening of reference ellipsoid 

A, C— equatorial and polar moments 
of inertia 

diurnal rotational speed 


No. 

Name 


U 

V 

W 

mm 

Tucson 

OPT 

-1996.0051 

-5042.6961 

3360.7748 


Pic du Midi 

OPT 

4686.1252 

11.6385 

4331.0499 

mm 

Johannesburg 

OPT 

5058.2628 

2698.0251 

-2799.8019 

— 

Goldstone 

RRT 

-2351.1949 

-4655.5944 

3661.0605 

— 

Woomera 

RRT 

-3978.5840 

3724.8986 

-3302.3278 

mm 

Johannesburg 

RRT 

5085.4787 

2668.3035 

-2768.7011 

■ 

Madri d 

RRT 

4845.7274 

- 360.0147 

4125.7615 


Moon 


m =81.3 


ratio of mass of the earth to that of the moon 


C2,0 

^ 2,2 


-.000207 

.0000207 


second-order harmonics 


= 3.6696 


ratio of equatorial semi -axis of earth ellipsoid to 
mean radius of the moon 


8 


.000629 - B = (C - A)/B ratio between moon's principal moments 

of inertia 
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“6 

Mascons In 10 Units of the Total Mass of the Moon 


Situated on the Surface of the Moon Sphere 


No. 

Mass 

Selenographic 

Longitude Latitude 

1 

+19.0694 

• -15 u 

15° 

2 

-20.0394 

20 

20 

3 

15.3805 

30 

- 5 

4 

-11.8910 

0 

-20 

5 

11.7378 

-45 

-40 

6 

-12.7844 

-50 

- 5 


2.8233 

-60 

30 


- 6.0431 

-15 

60 


10.0246 

50 

50 

10 

- 8.0542 

70 

10 

11 

8.1149 

60 

-30 

12 

- 8.3384 

30 

-60 


Trianqulation Points on the Surface of 


the Moon Sphere 
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Illuminating "Sun" 


4> n = 1.7399359 rad 

<f = .0172027913 rad/ day 

= <j> 0 + (T - T 0 ) 


constants in the equation giving the 
longitude of the illuminating vector 

T n - standard epoch 

1900.0 = 2415020.0 J.D. 

T - epoch in J.D. 
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3. POSITION DETERMINATION FROM EARTH-MOON OBSERVATIONS 
3.1 Lunar Far Side Positions from Apollo Trans-Earth 
Trajectory Photography 

1 2 
M.'D. Sprague and W. Riotte 

3.11 Sunmary of Theoretical Approach 

The investigation into the feasibility of extending selenodetic control 
on the far side of the lunar surface through an analytical photogrammetric 
method utilizing simulated Apollo trans-earth trajectory photography was 
completed [Sprague, 1971]. 

A brief summary of the theoretical approach to the problem, major 
conclusions and reconmendations are presented here. The theoretical 
procedure of the work was to generate synthetic photographic data of thirty 
triangulation points on the lunar surface far side as if taken with an 
80 mm Hasselblad hand-held camera from the Apollo command module during the 
initial phase of the vehicle's trans-earth trajectory. Such simulated 
"photos" were generated mathematically at one minute intervals over a 
forty-minute period (from TEI + 20 minutes to TEI + 1 hour) for four 
separate simulated Apollo missions. The data was generated from the 
simulated earth-moon environment program developed in [Papo, 1971]. This 

^Presently with the U.S. Army 

O 

Presently with the U.S. Air Force 



generated photo data, in the form of. photo coordinates of triangulation 
points and photo exterior orientation, plus survey coordinates of the 
triangulation points, was randomly perturbed within specified limits to 
simulate observational errors and then subjected to a block tri angulation 
adjustment calculation. Thus, adjusted values of triangulation point survey 
coordinates, photo coordinates and exterior orientation are obtained. 

Tri angulation points in the lunar equatorial band were used as known control 
points and their coordinates were constrained to ±300 meters in the adjust- 
ment. After adjustment, the adjusted values of the lunar survey coordinates 
for the tri angulation points appearing on photos used in the adjustment were 
compared with the corresponding true or theoretical coordinate to determine 
the degree of success of the method. 

The block tri angulation performed with theoretical data in this 
investigation differs from a conventional (earthbound) block tri angulation 
in two major respects. First, a conventional block normally has control 
points established around the perimeter of the block with some additional 
control within the block's interior. In this study, control points are 
established in a band or strip and appear through the center or along the 
edge of the block of photos depending upon the spacecraft trajectory. The 
second major difference lies in the pattern of photos in the block and the 
photo overlap resulting. For conventional blocks, the flight path and 
exposure timing is such that 60 percent end overlap and 30 to 60 percent 
side overlap is used for all photos within the block. Also, the photog- 
raphy is flown at approximately the same altitude above the terrain so 
that each photograph covers an area of nearly equal size. In the case of 
the lunar post-TEI photography, each successive photo is exposed with the 
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exposure station at a greater distance from the surface than for the 
preceding photograph. Therefore, each photo contains a larger area of the 
surface than the one taken before it. A point is quickly reached in the 
post-TEI trajectory when photo overlap is 100 percent (for one-minute photo 
interval), i.e., all points in the preceding photo are contained in the 
next photo exposed. 

3.12 Results of Simulated Block Triangulation Adjustments 

In the course of the investigation, block triangulation adjustment 
calculations were done with various sets of simulated photo data. The 
purpose was to experiment with numbers of photos used, camera "look angles" 
and time interval between "photos" used. Initially, various six-photo 
block adjustments were done for each of the four individual simulated 
trajectories. Then block adjustments were performed with increasing numbers 
of photo data sets up to a maximum of 16 photos both for single trajectories 
and for combinations of "photos" from two, three and four trajectories. 

Results of the triangulation adjustments in all cases indicated that 
some improvement was obtained in the adjusted location of a triangulation 
point compared to the initial perturbed or "unknown" location of that point. 
The degree of improvement was seen to be directly related to the number of 
photographic observations (number of photos the point appeared on), the 
number of trajectories or missions from which the point may be photographed, 
and the geometrical relationship between the point and the exposure stations 
The study revealed that the area of lunar surface on the far side for 
which this procedure may be effectively applied is severely restricted due 
to the very nature of the spacecraft trajectory after TEI. In order for the 
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vehicle to be placed in the proper trajectory to return to earth, the TEI 
burn must occur on the lunar far side near the 180° meridian, and by 20 
minutes after the injection burn, the spacecraft has proceeded in its 
trajectory to a point some 1200 to 1400 km from the surface and in the 
vicinity of the 100° to 120° E meridian. This, then, precludes the 
opportunity to photograph areas further east than approximately 130° to 135 E 
on the surface. Also, the position of the terminator may eliminate a large 
portion of area for photography depending upon the terminator location. 

The standard deviations for selenographic (lunar survey) coordinates 
of new extended control (NEC) points obtained in the results provides a 
means of comparing final adjusted value to the true or theoretical value of 
a coordinate. 

Inspection of the results of the block adjustments revealed certain 
similarities: (a) The path of the subvehicle point on the lunar sphere is 

the determining factor for the "observability" of a triangulation point 
(provided the point is illuminated), (b) A stronger solution (in terms of 
lower standard deviation and residuals) was obtained for those tri angulation 
points nearest the subvehicle path and thereby appearing on a maximum num- 
ber of photos included in the block adjustment. Points located at great 
distances from the track, such as the points at extreme northern latitudes 
and points east of the 120° E meridian, are not included in the photographic 
coverage until late in the period in which photography is taken. Thus, 
these points appear as fewer observations in the adjustment and the geometry 
of intersecting conjugate rays result in a weaker solution. 

(c) It was .noted throughout the study that standard deviations for 
Y-coordi nates of triangulation points was larger than for either the X or 
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Z coordinate- It was concluded that the reason for this was due to the 
location of the exposure stations with respect to the triangulation points. 

Fig. 3.1-1 illustrates the geometry involved. Primarily due to the curvature 
of the lunar surface, observations on points at substantial distances away 
from the subvehicle point introduce larger errors in the direction of the 
Y-coordinates than in either the X or Z direction for any small angular 
error <5. For points near the subvehicle point, such as point A in Fig. 

3.1-1 (a) and (b), a small angular error in the direction of the optical 
ray produces small errors SZ and 6X in position as shown. However, for 
points at increasingly larger distances from the subvehicle point as point 
B, the error in SY becomes dominant and is much qreater in size than <$X or 
sZ for points near the subvehicle point. This analysis holds true for 
exposures taken with the optical rays from the surface to the projection 
center lying primarily in the Y-coordinate direction as was the case for all 
four of the simulated missions. 

(d) In all cases as the number of photos in the block adjustment was 
increased, the values for standard deviations of survey coordinates 
decreased. Table 3.1-1 presents a comparison of the standard deviations for 
five representative NEC points resulting from four different block adjust- 
ments with simulated photo data from the Apollo 12 mission. 

In general, for sixteen-photo adjustments of photo data from single 
simulated missions, standard deviations for X and Z selenographic coordinates 
of NEC points were under 1.0 km, while the standard deviations of Y-coordinates 
were under 1.5 km. 

(e) Finally, results obtained through block adjustment of combined 
photo data taken from two or more simulated missions were much more promising 




Table 3.1-1 


Comparison of Standard Deviation for Five Representative NEC Points 
from Block Adjustments of the Simulated Apollo 14 Mission. 


NEC 
Pt. # 

Std. 

Dev. 

6 Block 
Block (1cm) 

10 Photo 
Block (lan) 

12 Photo 
Block (lan) 

16 Photo 
Block (km) 


ax 

0.550 

0.536 

0.497 

0.248 

12 

cry 

1.121 

1.088 

0.977 

0.563 


az 

0.218 

0.202 

0.184 

0.132 



(6) 

(10) 

(12) 

(16) 


ax 

1.427 

1.042 

0.925 

0» 844 

15 

ay 

2.335 

1.790 

1.549 

1.384 


az 

0.617 

0.497 

0.425 

0.378 



(5) 

(9) 

(11) 

(13) 


ax 

0.232 

0.154 

0.157 

0.122 

17 

cry 

1.132 . 

0.865 

0.750 

0.482 


az 

0.330 

0.264 

0.232 

0.173 



(5) 

(9) 

(11) 

(14) 


ax 

1.608 

0.326 

0.317 

0.192 

20 

ay 

4. 041 

1.139 

0.992 

0.603 


az 

1.253 

0.441 

0.378 

0.261 



(3) 

(7) 

(9) 

(12) 


ax 

3.169 

0.957 

0.833 

0.762 

23 

ay 

6.613 

2.522 

2.017 

1.738 


az 

2.209 

0.978 

0.761 

0.644 



(2) 

(5) 

(7) 

(9) 


(Numeral in parenthesis is the number of photos on which the point 
appears. ) 

than results secured through adjustment of single trajectory data. The 
sixteen-photo block adjustments performed using photo data from all four 
simulated missions achieved the strongest solution in that standard devia 
tions of coordinates were uniformly low for all NEC points appearing on 
photos from at least two of the four missions as shown in Table 3= 1-2. 
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Table 3.1-2 


Standard Deviation of Selenographic Coordinates (in meters) 
for NEC Points Resulting from Sixteen-Photo Adjustments 
of Combined Mission Data 


he c 

Point 

Ho. 

Ho. of 
Traj. 
Cont. 
Observ. 

Job 1245 
(TEE 25 
to 49) 
4 4min 

Job 1246 
(TEE 32 
to 56 ) 
4 4min 

Job 124? 
(TEE 22 
to 60 ) 

2 min * 

Job 1248 
(TEE 30 
to 60 ) 
2 min 

Max. 

Diff. 

(m) 




134m 

101 m 

137m 

36 

12 



237 

207 

241 

34 

<Tz. 

i * | 

asa 

113 

98 

113 

15 

c* 


HKH i 


B 1 1 1 11 

HBBSSIsI B 


13 cr Y 

4 

ImmsM 

■BB 

212 

216 

■ 

Oz. 


WKmm 


86 

101 




mM&gm 

612 

691 i 

669 

84 

15 

2 


799 

941 

908 

154 



Hhw 

252 

295 

277 

51 

HE3 


99 

121 

"^T 62 

122 

23 


4 

168 

226 

168 

228 

60 

BE! 


96 

108 

102 

108 

12 

B5I 



■nraa 


■K'» 

13 

17 <r«f 

4 


■ 

WE!m 


32 

(T* 




■SB 

MwiM 

15 




203 

wmm&rm ■ 

HU 

HH 

19 cv 

3 


394 



»F?jB 

fl*. 


Bu£9 

183 


BB 

■B 

<r* 




B B@!^B 

hush 

25 

20 <r Y 

4 

^EEIE 


247 


64 

<Ti 


MBIB 

BE1SB 

127 


14 

<r* 


13 a 

152 

■Htaa# 

■ 

6l 

23 

4 

30 6 

324 

KS: 


87 

<Tx 


130 

133 

■BOB 

IBB 

27 

|HB3 

■ 

HE3SH 



HH 

jHBi 


1 

B3S11 

B^lfc 



:^B 

BB 

1 

HEf&B 

BESB 

BB£B 

■BS 

BOB 

<TV 




■h£H£ 


1352 

28 «y 

1 

■ if j 1 


1 ho 


6372 



■K&S 

BBEB 

BESS 

B3££9 

1048 

■Ka 

ihkhk 


wm&sam 

—a«H 


KB 

II Eg 

2 


^BBs 

Bttt* 


|9 

BB 





BBS 

BIB 

ov 


302 

475 

406 

■ 


30 «v 

1 

1636 

4834 

1846 

B£kkl 


<T* 


385 

909 

441 

BmI 

BUBi: 
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This theoretical investigation has provided substantiation for the 
premise that the photogrammetric method of block triangulation through 
least-squares adjustment of post-TEI photos of the moon can provide a 
means of significantly improving selenodetic control on the lunar surface. 

The basic assumptions used throughout the procedure are listed below 
and they may also serve as limiting criteria for processing real data in 
preparation for block adjustment of actual mission photographs. 

(1) Standard deviation of photo coordinates (la value) is 10 micrometers. 

(2) Standard deviations for coordinates of those tri angulation points 
used as control within the block are 300 m or less. 

(3) Initial approximations for coordinates of NEC triangulation points 
under investigation are available with an accuracy of ±10 km. 

(4) Initial approximations for the elements of exterior orientation may 
be determined (by single-photo resection) to within 100 km for the 
exposure station coordinates and to within ±10° for the angular 
elements. 

With the above presuppositions in force, it is expected that improved 
selenographic coordinates for triangulation points may be obtained with 
resultant standard deviations of the order of 300 m to 1000 m if photo data 
from four or more Apollo missions with widely separated trajectories are 
subjected to the block triangulation adjustment. 

The conditions under which the optimum solution was achieved using 
simulated data were: 

(1) Photo data was extracted for use from four different simulated 

trajectories or missions which had widely separated subvehicle tracks. 
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(2) Four "photos" were selected from each of the four missions at four- 
minute intervals to comprise a sixteen-photo block. 

(3) Photos used were spread uniformly over a time period of thirty minutes 
from TEI + 20 minutes to TEI + 50 minutes. 

It is felt that the application of this method to post-TEI photography 
of the Apollo missions used in conjunction with the data reduction of 
other observational systems can contribute substantially in improving the 
existing selenodetic control net. 

3.13 Results from Apollo 15 Trans-Earth Photography 

The investigation towards the improvement of selenodetic control on 
the lunar limb through an analytical photogrammetric method utilizing 
Apollo 15 trans-earth photography was also completed [Riotte, 1972]. A 
brief summary of the problem, procedure, conclusions and recommendations 
are presented here. 

It was shov/n earlier that it is theoretically feasible to extend 
control from points of known location to unknown points on the lunar surface 
through the use of photogrammetric techniques with Apollo metric photography 
taken after the spacecraft leaves the moon's orbit. This project answers 
the question "Can this be accomplished with real data?" 

The successful Apollo 15 mission launched on July 26, 1971 was the 
first to carry a complete set of mapping cameras. The SIM (Scientific 
Instrument Module) contained among other items the Fairchild 3-inch mapping 
camera, a 3-inch stellar camera, and the ITEK 24-inch optical panoramic 
camera. This report used only the film from the 3-inch mapping camera. 



The metric camera's 4-1/2" x 4^1/2" format contained an array of 121 
reseaus and the calibrated coordinates of the reseaus were available from 
the calibration report. 

During the lunar orbits the metric camera took continuous pictures of 
the lighted portion of the lunar surface and continued taking pictures 
after the TEI (Trans-Earth Injection). It is this series of pictures with 
the moon receding from view that was used in this project. A formated 
picture is shown as an example in Fig. 3.1-2. There are several factors 
that should be noted from an evaluation of this sequence of pictures: 

(a) The trajectory during TEI is well below the equator; (b) there is a 
narrow band (approximately 20° longitude, 85°-105° E) of visible features 
bordered on one side by the highly reflective surface and on the other by 
the terminator; (c) the base height ratio between pairs of photographs is 
extremely small (approximately 0.04) and well below the desired limit of 
0.3 thus giving rear vertical photographs; (d) all measured photographs 
show 100% overlap; (e) the control points are found only on the photographed 
lunar limb; (f) the NEC (New Extended Control) points are located in the 
nadir region of the photographs. 

Identifying known control positions on the Apollo 15 film proved to 
be a time consuming task. The enhanced Orbiter IV photographs of the front 
side with the ACIC control points identified were examined and one control 
point was found. Three other control points were found from the pictures 
of the Lunar Landmark Control features. This is a most important network 
of control. It encircles the lunar equator and was established by reductions 
of numerous sightings of features by orbiting Apollo astronauts. Without 
this control network the project of extending control photogrammetrically 
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(p = - 31?. 29 X = 102°. 78 Range = 1915.1 (km) 



would be seriously hampered. 

Seven NEC points were selected as the points for which improved 
coordinates were desired. These are relatively small craters; however, 
they should be easily identifiable on pictures taken of the same area on 
succeeding missions. 

These eleven points on twelve photos were measured on a first-order 
instrument, an AP/C (Analytical Plotter/Commercial). The least squares 
reduction of 'the measurements to the photo system contained the computation 
for radial distortion. The radial distortion ranged from 0.8 urn to 44.4 ym 



and averaged 15.48 pm. This is within the range established in the 
camera specifications. 

These reduced observations of visible features, the calibrated focal 
length, the estimated values of the selenographic coordinates for each 
control and NEC parameter, the estimated values of the elements of exterior 
orientation parameters and the weights for the photo observations and 
estimated parameters formed the input to the FORTBLOCK adjustment program. 
The FORTBLOCK adjustment triangulation program performed a simultaneous 
least squares adjustment on the estimated parameters (i.e., elements of 
exterior orientation and lunar coordinates) based on the above inputs and 
the collinearity condition. 

The block adjustment was processed in 6 and 12 photo blocks. The 
program normally iterates internally three times; however, due to the 
unusual conditions of receding photography, the 6 and 12 photo blocks were 
iterated three times and then also for an additional three times. Table 
3.1- 3 shows the summary of the adjustment with the smallest standard 
deviations for the NEC points for the 12 photo block with six iterations. 
The adjusted values for the NEC points were then processed to provide the 
adjusted latitude, longitude and heights above or below a sphere of 
1738.1077 km radius. 

In summary (a) the standard deviations of the NEC points decrease as 
the number of photos and iterations are processed in the adjustment. The 
small differences in the standard deviations indicate the solution has 
reached its limitation with the six iterations of the 12 photo block. 

(b) The standard deviation for the Y coordinate for NEC points 12-17 
which are located near the nadir region of .the photograph are always 
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Table 3.1-3 


Summary of Results 12 Photo Block Adjustment {6 Iterations) 


NEC 
Pn i nf- 

No. of 
Photos 

Adjusted Values of 
Survey Coordinates (km-) 

Residuals (km) 

Standard Deviation (km) 

No. 

App. On 

x a 

Y a 

Za 

v* 


Zx 

a. 


CTx 

11 

12 

- 100. 378 

172.3, „743 

-• 133.725 

- 1. 879 

7. 369 

- ' . 121 

.538 

. 810 

. 928 

12 

8 

- 152. 632 

1438.456 

- 950.355 

- 6.400 

3.450 

- 9.471 

1. 291 

2.529 

1.558 

13 

12 

- 119, 109 

130G. 499 

- 1132.546 

- 5. 609 

- 11.190 

- 15. 110 

1. 389 

2.922 

1.757 

14 

12 

123. 220 

1327, 734 

- 1104.417 

2.454 

4.409 

- 6.331 

1.494 

2. 694 

1.526 

15 

12 

154, 426 

1083.855 

- 1347.756 

- 4.039 

4.989 

.320 

' 1. 611 

3. 602 

2.223 

16 

8 

- 1. 244 

.-970. 962 

- 1439. 079 

- 8.006 

- 3. 084 

- 5.761 

1.632 

4.253 

2, 658 

17 

' 6 

- 224. 673 

1373.941 

- 1028, 307 

2. 055 

- 1,446 

- 1,330 

1.311 

2. 852 

* 

1.753 


a o “ 1. 56 

Average p xy “ - . 07504 
Average p„ = - . 14994 

Average p y2 = 


57366 


























greater than the standard deviation of X or Z. The correlation coefficients 
p and p were greater than p and also these NEC points are relatively 

Aj j A AA 

highly correlated with each other. This can be attributed to the poor 
geometry of intersecting rays to points near the center of the photograph 
and also because all the photographs were taken with the selenographic 
Y axis nominally toward the spacecraft, thus the convergence is not as 
precise as in the X-Z direction. 

(c) The adjustment of the 6 and 12 photo block had an effect on the 
residuals and standard deviations of the four control points. In almost 
all cases the resulting standard deviations increased slightly from the 
values provided by NASA. It is felt that the unique geometry not only 
from the spacecraft traversing away from the lunar surface but also the 
location of the control in just the northern limb of the photographed moon 
creating a very narrow cone of intersecting rays contributed to this 
problem. 

3.14 Conclusions and Recommendations 

The results show that the system originally specified in theory 
[Sprague, 1971] is workable with real data. The location of relatively 
unknown features on the lunar limb can be improved from 20 km to 1-3 km 
by use of photogrammetric procedures with metric film taken on one Apollo 
trans-earth trajectory. As a result of the investigation (a) it was 
reconwended that the same project be accomplished using the data and film 
from Apollo 15, 16 and 17. The solution using photographs of the same area 
from different post-TEI trajectories will provide improved solutions 
through additional perspective rays. 
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(b) The TEI trajectories should be aligned in a more equatorial region. 
This would allow all the control points within the region to be observed 

and measured and would allow all the control points to appear in the central 
region of the photographed moon. 

(c) The Lunar Landmark Control network should be densified with 
additional sightings on succeeding Apollo missions. Improved solutions 
of extension of control problems using photogrammetry can be provided 
when additional control is available. 
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3.2 Lunar Near-Side Positions from Earth-Based Photography 


H.B. Papo^ 

3.21 Introduction 

The problem to be treated in this paper can be defined broadly as 
mapping the moon . By mapping the moon we mean here the determination of 
positions of specific points and features on its surface relative to a 
certain coordinate system. Coordinate systems can be chosen at will, but 
there is added importance to those which are physically significant. For 
the moon a Cartesian coordinate system can be defined which is centered at 
its mass center and is oriented along its three principal axes of inertia. 
The above system will be referred to as the selenodetic coordinate system 
of the moon. 

There are in existence many different maps of the moon which define 
coordinates of features on its surface. Unfortunately, the same features 
have different coordinates on different maps. The main reason for these 
differences is the lack of a unique fundamental control network on the moon. 
A network of 20 or 30 distinct features on the visible side of the moon 
which have precise coordinates with respect to the selenodetic coordinate 
system could serve as fundamental control on the moon. In order to avoid 
the large discrepancies between the various maps, any subsequent densifi ca- 
tion of control for mapping purposes should be derived from the same 
fundamental control points. Thus the specific problem to be treated in 

• i 

Presently at the Technion, Haifa, Israel. 
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this paper is the establishment of a fundamental selenodetic control network 
on the moon. 

A basic difference between the earth and the moon in establishing 
such a control system is in the fact that there are no observations of 
geodetic significance which have been performed from the lunar surface. 

On earth we use classical or satellite triangulation and trilateration 
geodetic astronomy and gravimetry in all of which the observations are 
conducted from points on the earth's surface. On the moon such an analogy 
is still nonexistent. The only observational material at our disposal are 
remote observations taken directly from the earth's surface or indirectly 
from a spacecraft which is tracked continuously from earth. Because of 
this important difference, establishing control on the moon is inseparable 
from two aspects of its motion in space, namely, the instantaneous position 
of the center of mass of the moon with respect to the earth and the 
instantaneous orientation of the selenodetic coordinate system with respect 
to inertial space. These two aspects are by themselves rather difficult 
problems in dynamical astronomy; however, there is no way of establishing 
a control of physical significance without insuring its consistency with the 
motion of the moon. 

The numerical lunar ephemeris introduced by the Jet Propulsion Laboratory 
under the code name LE-16 represents the translatory motion of the lunar mass 
center with respect to the geocentric 1950 mean equatorial coordinate system. 
According to the definition given above, the instantaneous position of the 
origin of the selenodetic coordinate system at a specific epoch is identified 
in our study with the LE-16 lunar ephemeris. 

The second aspect of the motion of the moon, namely, its rotation 
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about an axis passing through its mass center, was made an integral part 
of our solution for fundamental control. At any specific epoch the rotation 
of the moon defines the instantaneous orientation of the selenodetic 
coordinate system with respect to inertial space. In our study the parameters 
of orientation of the selenodetic system are identified with the physical 
libration of the moon. A full account on the approach taken in solving the 
physical li brat ions of the moon is given in section 2.2. 

The solution for fundamental control presented in this paper is based 
on remote optical observations of the lunar surface. Photography of the 
moon taken from the earth or from satellite-borne cameras as well as direct 
angular observations to specific features on the moon taken from the earth 
or from a satellite all form the main bulk of observational data. Auxiliary 
data designed to provide orientation and scale to the optical data is also 
considered (stellar camera, star background lunar photography, etc.). 

Another important although "nonoptical" auxiliary data source are 
range-rate observations from the earth to a spacecraft involved in optical 
observations of the moon. As the instantaneous position of the spacecraft 
with respect to the moon at any specific epoch is derived from analysis of 
range-rate observations, this data is also processed in the general solution 
for fundamental control . 

3.22 Basic Concepts 

For a fundamental orientation frame of reference, a hypothetical 
inertial coordinate system is considered which coincides with the ecliptic 
mean coordinate system at some arbitrary standard epoch and is defined as 
a Newtonian Frame of Reference [Brouwer and Clemence, 1961, p. 3]. The 
orientation of any other Cartesian coordinate system with respect to the 
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inertial system is defined through three Eulerian angles necessary to 
rotate the particular Cartesian system into the inertial or vice versa 
(see Fig. 3.2-1). The term inertial coordinate system, to be denoted by 
XYZ, is used for any coordinate system (having an arbitrary origin) which 
is parallel to the fundamental orientation frame. Thus, one may have a 
geocentric inertial system, an inertial system centered at a satellite, 
etc. 

The orientation of a physical body is defined through a Cartesian 
coordinate system fixed to the body. Thus, for example, the orientation of 
the earth is defined through the so-called average terrestrial coordinate 
system which is fixed to the earth's crust and related to a geocentric 
inertial system by means of three Eulerian angles. 

There are three Cartesian systems which are of primary importance: 

(a) The average terrestrial system is centered at the mass center of 
the earth and is oriented with respect to the crust through the CIO pole 
and the so-called mean observatory [Mueller, 1969]. It is denoted by UVW. 

(b) The selenodetic system is centered at the mass center of the 

and is oriented along its principal axes. It is denoted by xyz. 

(c) The optical observations reference system is centered at the 
projection center (from which the optical rays emanate) and is generally 
oriented so that the primary axis (Bj) points towards the moon. It is 
denoted by B B B . 

I u O 

The general rotational motion of a physical body with respect to the 
inertial system is described by the change with time in the Eulerian orienta 
tion angles of the appropriate coordinate system (the one fixed to the body) 
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Fig. 3.2-1 Earth-moon environment. 


xyz 


V,\,€ - 
<P> M - 

o P 

o D 

T 


the geocenter 
the selenocenter 
a satellite 

inertial coordinate systems 

average terrestrial coordinate system 

selenodetic coordinate system 

reference frames for optical observations 

Eulerian orientation angles of the Earth 

Eulerian orientation, angles of the Moon 

observing station for optical observations 

tracking station for range and range rate observations 

triangulation point on the Moon 


range (range rate) between O 0 and S 


optical ray (direction) between O p or S and T. 




Thus, for example, in order to study the rotational motion of the moon one 
should consider the changes in the Eulerian angles (<j>,4>,9) which relate the 
oon-fixed (xyz) system to the inertial (XYZ) system. From this definition 
of rotational motion it is clear how rotation and orientation of a body are 
related to one another. As the body rotates, its orientation changes and 
by considering the instantaneous set of Eulerian angles at a particular 
epoch one has the means for defining the orientation of the body at that 
epoch. 

The rotation of a physical body is governed by a set of equations of 
rotational motion (second-order differential equations) in which the 
external forces acting on the body are represented. The solution of these 
equations of motion results in the Eulerian orientation angles as functions 
of time. 

The position with respect to the UVW system of points on the surface 
of the earth is defined similarly with respect to the selenodetic (xyz) 
system. 

The position and velocity (state vector) of points in space with 
respect to a particular inertial coordinate system (selenocentric, geocentric, 
etc.) are defined by Cartesian coordinates (XYZ). Thus, the geocentric 
state vector of the moon's mass center is given by the LE-16 ephemeris 
[O'Handley et al., 1969] where the inertial XYZ system has been defined as 
identical to the mean equatorial system of 1950.0. 

All optical observations of the moon (photographs, direction measure- 
ments, etc.) are treated as light rays emanating from a projection center 
whose selenocentric inertial coordinates are known or are being estimated 
in a least squares process. The individual ray from a bundle (the rays 
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emanating from the same projection center form a bundle) is related to a 
reference optical frame B^B^ by two angular quantities. As mentioned 
above, the is related to the XYZ system centered at the projection 

center by three Eulerian angles or in general by an orthogonal transformation 
matrix. 

If the projection center is on board a spacecraft, range and range- 
rate data from the earth to the spacecraft serve for determining its 
selenocentric orbit. The ordinary orbit determination procedure solves 
for the state vector (position and velocity) at a standard (initial) epoch 
as well as for a number of constants. As some of these constants are 
dominant factors in the rotational motion of the moon and also as the 
projection center from which the optical observations are made lies along 
the trajectory of the spacecraft, it is necessary to process the range 
and range-rate tracking data together with the optical data and thus obtain 
a solution for the trajectory and for the relevant constants which is 
consistent with both types of data. 

To summarize, the moon is regarded as a rigid celestial body rotating 
in space under the gravitational influence of the earth and the sun while 
around it there is a variety of sensors engaged in optical observations. 

The optical observations as well as the range-rate and other auxiliary data 
are modeled into observation equations in terms of the parameters of the 
solution and are treated thereafter by a weighted least squares procedure. 

3.23 Adjustment Theory 

The mathematical model for an optical observation is derived from the 
triangle formed by the selenocenter M, the particular control point being 
observed T and the projection center S or Op (see Fig. 3.2-1). 
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If we denote the vectors connecting the three vertices M,T,S according 
to Fig. 3.2-2, the mathematical model simply is 

F = u + s- t “ 0 

where 

u represents the optical observation 

t,s represent the selenocentric position of the control point and the 
projection center respectively. 

The mathematical model for range-rate observations is derived from the 
triangle defined by M, as before, the spacecraft S and the tracking station 
on earth 0 (see Fig. 3.2-3). 



If R denotes the instantaneous position of S with respect to 0 and R is the 
corresponding time derivative of R, the following equation holds: 

H = -jij (R - R) - H b - o 

where r is the range (distance) US and r is the observed rate of change of r. 
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The above nonlinear models have to be linearized before an adjustment 
procedure can be developed. By linearization we mean the usual Taylor 
expansion of the models in terms of small corrections to the approximate 
values of the parameters. The parameters considered in the adjustment are 
of two basic types, namely, parameters which are being sought to be referred 
to as permanent and parameters which are not of interest but which are 
indispensable in modeling the observations. The permanent parameters of the 
solution are placed into three groups: 

(1) Cartesian coordinates of the control points (3 x n) 

(2) parameters of orientation of the moon (6) 

(3) selected physical* parameters of the moon (3). 

If we choose to limit the number of control points to say 22, the total 
number of unknowns would be 75 which is not too large for a simultaneous 
solution on a medium size computer. The situation is enhanced further by 
the fact that the normal equations can be created by the summation of 
individual layers. The final set of the 75 normal equations is obtained by 
sequential addition of the contribution of numerous batches of optical 
observations. In this way the maximum size of the computer program required 
to form and solve the normal equations is determined mainly by the number 
of control points and not be the number of observations which in our case 
may run into the thousands. 

As the physical librations form an integral part in the solution, the 
particular logic of the adjustment procedure was developed with the numerical 
integration program for the physical librations serving as a backbone. The 
sequence of steps is as follows: 

(I) The integrating program is operated up to the "next" epoch at which 
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a bundle was observed. 

(2) The integration is arrested and the information generated by the 
integrating program for that epoch (physical libration angles, state 
transition and parameter sensitivity matrices) together with informa- 
tion gathered by the optical bundle are used in evaluating the partial 
derivatives as developed in [Papo, 1971, section 2.4]. ■ 

(3) A layer of the normal matrix, the constant vector and the considered 
parameters contribution matrix is generated and added to the corres- 
ponding matrices where layers from previous bundles have been 
accumulated. 

(4) The integration is resumed until the epoch of the "next" bundle and so 
on, until all the bundles in the batch have been processed. It should 
be remembered that the normal matrix generated is that of the permanent 
parameters only, the solution for the auxiliary parameters being 
"folded in." 

(5) The a priori covariances of the permanent parameters are added to the 
normal matrix followed by inversion and subsequent evaluation of the 
solution vector (corrections to the starting values of the permanent 
parameters) and the full covariance matrix of the corrected (adjusted) 
parameters. 

The above procedure is used for processing earth-based optical observations. 
Processing satellite-based optical observations is done in a similar way only 
now the range-rate data is added in the process and also the number of 
auxiliary (nonpermanent) parameters is increased to allow for the orbit 
determination of the satellite. 


70 



3.24 Numerical Experiments 

In order to test numerically the mathematical procedure developed in 
this study, a simulated environment was created which reflects very 
closely the true world. The earth, the moon and a variety of satellites 
move and rotate in this simulated environment strictly according to the laws 
of Newton and Kepler. 

The solution as proposed was tested through the analysis of simulated 
optical observations. The observational material was generated free of any 
unaccounted phenomena and simulated very closely real observations. 

Thirty well distributed points on the front side of the moon were 
chosen as the control network and over a period of one year thirty simulated 
photographs (bundles of optical rays) were taken from observatories on 
earth (see Fig. 3.2-4). Solutions were performed varying-the accuracy of 
the observational material and also varying the initial (approximate) 
values of the parameters. 



Fig. 3.2-4 Fundamental control network. 
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The objectives of the experiment were twofold: 

(1) To determine the quality of control solution which can be obtained from 
the simulated observational material. This was done by inspecting the 
covariance matrix of the solution vector of the parameters. 

(2) To find to what extent the adjustment procedure is capable of recover- 
ing shifts which were introduced in the nominal values of the parameters. 
This was accomplished by using the shifted values as first approxima- 
tions in the adjustment and by subsequent inspection of the solution 
vector. 

Table 3.2-1 presents diagonal submatrices of the covariance matrix of 
the solution for control coordinates in kilometers squared. It also gives 
the corresponding correlations. A clear pattern of higher correlations can 
be observed as follows: the x coordinates are correlated among themselves 
and the same applies to the y and z coordinates respectively. A possible 
reason for this phenomena may be that the orientation of the moon was made 
a part of the solution. Typical for all the points is that the x coordinates 
are determined much less accurately as compared to the y and z coordinates. 
The main reason for the apparent deficiency in determining the x coordinates 
is undoubtedly the poor geometry of any optical observation of' the moon as 
taken from the earth. We should realize that the maximum angle of conver- 
gence between optical rays from different bundles is usually much less 
than 23 degrees (see Fig. 3.2-5). The average standard deviations for x, 
y and z are 700, 150 and 150 m. 

Table 3. 2-2 presents the correct coordinates (nominal) of the control 
points and the negative values of the solution vector after one iteration 
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Table 3.2-1 


Covariance and Correlation Matrices for Solution of Network I 
from 30 Bundles. Triangulation Points. 



Parts of Covariance/ Correlation Matrix on this Page 


11 12 13 

0.483 0.C09 — 0.001- 0.303-- 0.QG9- -0.000 — 0.2-71- -0 .0C1 -0-.001 r 

O. C61 | 0.047 C.000 -0.019 0.046 0.000 -0.C37 0.C18 0.0C1 11 

— 0. 008— O.OO-i-j— 0-.0-1-5 — &.Q01 0.000—0.012 - 0.001— -O.COC - 0.010 

0.624 -0.128 0.016 1 0.488 -C.G23 -0.001 0.309 -0.007 -0.003 . 

0.060- - 0.953 -- C.0C-2~-8. I46 { - C.050- - C.OOO— Q.C39- C. C-19---0 .00 1 12 

-0.004 0.018 0.818 -0.009 0.018 1 0.0 14 -0.001 Q.000 0.011 

-- 0.529 — 0.22-5 — C. GC-8 — 0 . 600 — -0. 2-3 8--0.012 j Q . 544- —0 .012 -0. 002 
-0.009 0.679 0.025-0.077 C.669 0,015-0.135 1 0.015 -0.000 13 

G-*rOO-8 — 0-. 02-3 — 0» 664 — — 0-. 03-1 — 0.022 — 0-« 7 5 5 — 0. 023---— 0-. 0-10-j — 0 .0-1-5- - 
xyzxy zxyz 

, _ . jg ~ 19 

- C . 48 0 -0 . 0€-7 — € . OG-1 — (T. 2 83 — 0.0 02 — 0.0 02 - 0.30 9 — 0 .00-7 — O-.O 0 1- 

-0.044 1 0.058 -C.OOO 0.040 0.023 -O.OC1 -G.G27 0.G50 O.COO 17 

— 0". 01r4— — 0-. 00-2- — G-» 0-1-4 — — 0-.00-1 — —0.-000— — 0- .0-12 O-.O 02 — — 0-. 000 0 . 0-1-2- 

C. 551 0.226-0.013 1 C.551 C.CC8 -0.005 0.261 0.035 0.000 

—0.022- -0 . 74G—-0 .0-13-— 0.083-1 — 0 .0-17— -0.000-— 0.-CC 8-- G.G22 — 0.000- - 18 
0. 028 -0.033 0.818 -C. 059 -0.030 1 0.016 0.004 -0.001 0.011 

- -- Q.6-39 — 0 . 1 6 2 — 0 . 024 0 .505 — 0.0 8-7-- 0.041 ] 0.485- - -0.027 — 0.000 

-0.045 0.954 -0.012 0.214 0.765 -0.045 -0.177 1 0.047 0.000 19 

— 0. GO-8 0.0-P7 — 0.-822 — 0.00-5 — 0.006 — 0.-780—^0.006— 0-.0-I4-]-0-.0 14->-- 


x y z X y z x y-z 

--- ^ ' ‘ ‘ .29 


0*490 — 0-.013 — G.OC-G — 0.300— 0. 019— 0. 001 — 0. 294- -0.010 •• 0.002 
0.144) 0.017 0.001 -0.005 0.018 -0.000 -0.000 Q.Q15 -0.000 27 

0.CC5- 0.0421 - C.014--0.000— 0.001- -0. 0-1-1 — 0.000 — 0*001- 0.011 
0,629 -0.059 -0.0C3 | 0.466 -0.007 0.002 0.257 -C.CC4 0.002 
0.161 -0.840 — 0.046--0.0621 - 0.028 --C .000— - 0 .GC4-- 0.018 -O.COl 28 
0.C17 -0.020 C. 800 0.025-0.017 ) 0.014 0.001 -0.000 0.011 
0.604- -0.003 C.0G6 0 .625 0.C34 -0. 014 | - 0.484- -O.COl 0.G01 i 

0.IC6 0.859 0.049 -0.048 0.8L3 -0 .026 -C .OlO l 0.017 -0.001 29 

0.021 -0.0 28 - C. 75 5 - G.02Q ~-0. C2 6 - -G.82 8-- 0.009 -0 ,033j 0.014 
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Table 3.2-2 

Solution Vector for Test (a) in Experiment (ii) 


Solution Vector 
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4 * 
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00 
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Fig. 3.2-5 From the earth the moon is seen through a 
narrow window. 

of the adjustment. The shifts in the coordinates were obtained by simply 
removing the decimal fraction of the coordinates. So the degree of recovery 
of the shifts is within 10-20 m in y and z which is remarkable considering 
the low quality of the solution as exhibited by the covariances. As could 
be expected after only one iteration, the recovery in x is much less 
efficient. 

3.25 Summary 

The main characteristics of the solution for selenodetic control as 
developed in this study are summarized in the following: 

(a) The solution is consistent with the motion of the moon in space. 

(b) The solution for the orientation of the moon in space, is part of 
the general solution. 

(c) Optical data obtained from the earth or from a spacecraft are 
processed uniformly, thus avoiding inconsistencies between solutions based 
on either of the two sources of data. 

(d) All the observations needed for the solution (optical, range, and 
range-rate) are processed simultaneously in a weighted least squares 
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procedure where the parameters are constrained according to their a priori 
covariances. 

(e) The adjustment procedure can be programmed for used with available 
electronic computers where the core size required and the computer time 
for processing the data are reasonable and make the application of the 
solution to processing real data a feasible proposition. 
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3.3 Positioning from Lunar Laser Ranging 


F. 


Fajemirokun 


and F. Hotter 


3.31 Introduction 

Lunar laser ranging has provided another means by which geodetic 
control can be established on the moon. It is now estimated that distances 
between stations on the earth and lunar retroref lectors can be measured 
with a precision of ±15 cm with laser techniques [Alley et a!., 1970]. • 

An improvement in this precision is expected in the future. 

As opposed to the methods of establishing geodetic control for earth 
mapping which involve measurements on the earth's surface itself, past 
selenodetic controls have been obtained mainly through the use of earth- 
based or lunar satellite-based observations. Thus, the determination of 
coordinates of points on the moon through these methods is rigorously tied 
to the following extraneous parameters: 

(1) the coordinates of the moon's center of mass in a geocentric inertial 
coordinate system, 

(2) the orientation parameters of an earth-fixed coordinate system with 
respect to a moon-fixed coordinate system, 

(3) parameters of the atmospheric refraction model used when reducing the 
earth-based observations. 

Lunar laser ranges of such precision, as mentioned above, expected not only 
to improve the accuracy of positions of points on the moon, but also to 

^ Presently at the University of Lagos, Lagos, Nigeria 
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contribute to more accurate determination of the geocentric positions of 
laser stations on the earth. Furthermore, it is expected that the parameters 
of the earth's and moon's orientation in space can be also improved upon. 

The relative position of a lunar point (such as the retroreflector) 
with respect to a station on the earth depends on a number of physical and 
geometric parameters of the earth, moon and the earth-moon dynamic system. 
Therefore, a measured distance between an earth station and a lunar point 
can be expressed mathematically as a function of the following parameters: 

a. the selenodetic coordinates of the lunar point, 

b. the geocentric coordinates of the center of the moon, 

c. the geodetic coordinates of the earth station, 

d. parameters of the orientation of a moon-fixed coordinate system 
with respect to a celestial coordinate system, 

e. parameters defining the orientation of the "average" terrestrial 
(earth-fixed) coordinate system with respect to the celestial 
coordinate system. 

Before deriving the earth-moon distance equations, the above-mentioned 
coordinate systems and orientation parameters will be given brief 
descriptions: 

The selenodetic coordinates — longitude i, latitude b, and radius 1 — 
define conventionally the location of a point on the lunar surface. The 
selenodetic latitude is measured from the lunar equator positive toward 
the north lunar pole. The prime meridian (z = 0) is a plane which contains 
the lunar axis of rotation and the earth-moon line at zero geometric 
libration. The longitude is measured from this prime meridian positive 
in the direction of rotation. The Cartesian selenodetic coordinate system 


78 



(x,y, z) is fixed with respect to the solid moon, is centered at the moon's 
center of mass (selenocenter) and is related to the coordinates z, b, and 
r the usual way [Mueller, 1969]. 

The geocentric coordinates of the selenocenter is obtained from a 
lunar ephemeris based on a particular lunar theory which attempts to -solve 
the differential equations of motion of a perturbed two-body system. These 
differential equations can be integrated either analytically [Brown, 1908] 
or numerically [O’ Handley et al., 1969]. 

The "average" terrestrial coordinate system (u,v,w) in which coordinates 
of the earth stations are expressed is defined by the Conventional Inter-, 
national Origin (w-axisj and the Greenwich Mean Astronomical Meridian as 
determined by the Bureau International de VHeure (u-axis) [Mueller, 1969]. 
The coordinate system is fixed to the solid earth and is centered at its 
center of mass. The Cartesian coordinates (u,v,w) can be computed from 
the conventional geodetic latitude (a), longitude (x) and height (h) the 
customary way [Mueller, 1969]. 

3.32 Numerical Experiments 

Numerical experiments were performed for the purpose' of investigating 
the expected accuracies of the parameters in the adjustment model using, 
laser observations of currently expected precision. Since no real data 
were available, laser distances were simulated for the purpose of these 
experiments. 

The numerical experiment utilized the three lunar retro ref lectors so 
far deposited on the moon in the Apollo program, to which distances were 
simulated from certain existing laser stations on the earth. First, a 
set of experiments were performed using only the laser station at the 



McDonald Observatory ranging to the three retroreflectors. Later, in order 
to find out the effect of multiple observing stations on the recovery of 
the parameters, two stations were added: Mt. Stromlo and the Crimean 

Observatories. Tables 3.3-1 and 3.3-2 contain the coordinates of the 
laser stations and the lunar retroref 1 ectors . 

For the limited objectives of these experiments, only the following 
parameters were included in the adjustment as unknowns: 

(1) The geocentric (u,v,w) coordinates of the laser stations on the 
earth. 

(2) The six parameters of the orientation of the (u,v,w) "average"' 
terrestrial coordinate system with respect to the mean ecliptic, 
coordinate' system of 1969.0. These are the -Eulerian angles and 
their time derivatives .for the epoch of 1969.0. 

(3} The selenodetic (x,y,z) coordinates of the lunar retroreflectors. 

(4) The parameters of the orientation of the {x,y, 2 ) lunar coordinate 
system with respect to the mean ecliptic system of 1969.0, repre- 
sented by the physical libration parameters— and their 
derivatives at the epoch of 1969.0. 

(5) Three physical parameters of the moon given by the spherical 
harmonic coefficients C £Q , C 22 and the moment of inertia ratio 



The simulated data were assumed to be completely free of systematic 
errors, and the geocentric coordinates of the moon's center of mass were 
assumed to be correct as given in the ephemeris [O'Handley et al., 1969]. 
The numerical integration of the earth's Eulerian angles is outlined 
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Table 3.3-1 


Geocentric Coordinates of Laser Stations 


Point No. 

u(lon) 

v (km) 

w(kni) 

Location 

1 

-1330. 81462 

-5328.78935 

3235.69752 

Ft. Davis, Texas 

2 

-4466.54586 

2683.24104 

-3667.44266 

Mt. Stromlo, Australia 

3 

3784.28692 

2552.21344 

4440.46175 

Crimean, USSR 


Table 3.3-2 

Selenodetic Coordinates of Lunar Retroref lectors 


Point No. 

x(km) 

y(km> 

z(km) 

Remarks 

1 

1591.42945 

691.97322 

19.26380 

Apollo 11 Mission 

2 

1652. 86823 

-520.44969 

-110.54128 

Apollo 14 Mission 

3 

1554.85951 

99.26753 

762.35341 

Apollo 15 Mission 
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in [Fajemirokun, 1971]- In the Encke-type of integration, the quantities 
integrated are the following: 
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tion matrix U, which was obtained by numerical differentiation is given as 
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For the corresponding lunar parameters, the method described in 
[Papo, 1971] was used which integrates the physical libration parameters 
(t,o,p ,t,ct and p). The Eulerian angles and their time rates were obtained 
from the physical libration quantifies and the Cassini expressions as follows 
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where 


L is the mean longitude of the moon in its orbit 
q is the longitude of the mean ascending node of the moon's orbit 

I is the mean inclination of the lunar equator with respect to the 

ecliptic 

In addition to the physical libration quantities and the state transition 
matrix, the integration program also computes the parameter sensitivity 
matrix 

3[t a p i a p] 

5 _ — - 

a[c g c y 

22 20 

and the state transition matrix 

r . . .,T 

8 It o p t <j p I 
U* = J: » 

The adjustment method regards all the parameters (in addition to the 
simulated distances) as observations with associated weights. The relative 
weights of the parameters were computed as the inverse of their estimated 
variances, thereby choosing the variance of unit weight to be equal to one. 

The values of the variances of parameters were selected to conform 
with the level of uncertainties in the present knowledge of these quantities. 
The chosen standard deviations of the parameters are given in the last 
column of Table 3.3-3. 

3.33 Results 

The purpose of this experiment was to investigate the internal 
precision of the adjustment system. There were three cases investigated: 
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(1) Laser distances were simulated between McDonald Observatory and the 

three lunar retrorefl ectors . The number and period of observations 
chosen were as follows: 50 observations over a three-month period, 

100 observations over a six-month period, 150 observations over a 
nine-month period, and 200 observations over one year. 

(2) The number of observations was held to 50 and 100 while the period of 
observations was varied. 

(3) The simulated distances used for the adjustments comprised of distances 
between all the three laser stations and retrorefl ectors. Number and 
period of observations chosen were the same as in (1). 

Table 3.3-3 gives the diagonal elements of the variance-covariance matrix 
for case (1). The general pattern shows, as can be expected, that the param- 
eter variances decrease as more observations are made over longer periods of 
time. It can also be seen that in addition to a general improvement in lunar 
position, we can also expect an improvement in the determination of geocentric 
positions of stations on the earth's surface. A general improvement in the 
accuracy of determining the orientation parameters of both the earth and the 
moon, and the moon's lower order spherical harmonic coefficients can also be 
expected. Exceptions seem to be the variances of the physical libration in 
the node and its time derivative (a and a). Even with 200 observations over 
a one-year period, the variance obtained for a and a are 250'.'7 and 29.5 sec / 
day 2 , respectively, compared to the a priori values of 400" and 100 sec /day . 
A good explanation for this poor determination of a and a was not found, 
especially since the dominant term of libration in node has a period of about 
one month. On the other hand, it should be borne in mind that the node it- 
self is rather poorly defined through the intersection of two planes at an 
angle of only 1 ?5. 
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Table 3.3-3 


Variances of Parameters for Case 1, (One Laser Station Observing, 
Number and Period of Observations Varied) 



50 Observations 

100 Observations 

150 Observations 

200 Observations 

A Priori 

Parameters 

3 Months 

6 Months 

9 Months 

12 Months 

Variance 


meter 2 

meter 3 

meter 3 

meter 3 

metei' 3 

u 

312.7 

1.0 

0.3 

0.08 

625 

V 

19.6 

0.07 

0.02 

0.01 

625 

w 

36.0 

8.9 

1.8 

1.4 

625 

*1 

18.3 

2.9 

1.2 

0.5 

10 s 

y x 

57.1 

10.1 

4.2 . 

1.9 

10 6 

Z 1 

646.3 

220.9 

• 85.8 

57.7 

- 10 s 

X 2 

10.5 

2.3 

.1.1 

0.4 

' 10 s 

y 3 

138.7 

27.6 

10.6 

4.4 

10 6 


627.4 

191.9 

66. 6 

50.1 

10 s 


17.5 

14.4 

10.1 

6.0 

10 6 

y 3 

113.2 

26.7 

9.8 

6.3 

10 s 

z s 

539.6 

171.5 

63.2 

44.6 

10 s 


HQSSHHH 






0.6 X 10“ 3 

0 . 1 X 10" 3 

0.2 X 107 4 

0.1X10" 4 

1.0 

6 2 

0.3 X 10" 3 

0. 2 X 10~ 3 

0. 8 X 10" 4 

0.6 XlO" 4 

1.0 

63 

0.5 

0. 1 X 10" 2 

0.3 XlO -3 

0.1X 10" 3 

1.0 

r 

1.5 

0.3 

0.1 

0.04 

400 

a 

397.8 

382.0 

318.7 

250.7 

400 

p 

1.4 

1.2 

0.8 

0.4 

400 

C 22 

0. 3 X 10" 2 

0.2 xlO" 3 

0.5 XlO" 4 

0.2 xlO" 4 

0.25 

JS 

0.4 X 10" 3 

0.3 xl0~ 4 

0.8 Xl0" s 

0.3 XlO" 5 

4.0 

C® 

0.1 X 10" 3 

0. 1 x 10~ 3 

■0. 1 X 10~ 3 

0. 1 x 10" 3 

0.1x10* 




mu yiiv 

mmTTJmMm 

lscc 2 /dav 

Si 

0. 2 x 10" 2 

0. 2 X 10*- 3 

0.1 X 10** 3 

0. 8 X 10" 4 

0.25 

6s 

0. 6 X 10" 1 

0. 7 X 10“ 2 

0.6 X 10" 3 

0.3 XlO" 3 

0.25 

6 6 

0.4 X 10" 3 

0.5 XlO" 4 

0.4 X 10" 5 

0.2 X10" S 

0.25 

* 

r 

0. 8 X 10“ 3 

0. 6 X 10" 4 

0. 1 x 10" 4 

0.6 XlO" 5 

100 

CT 

98.2 

81.3 

55.1 

29.5 

100 

• 

P 

0.02 

0.015 

0.011 

0.009 

100 
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Another trend that could be noticed in Table 3.3-3 is that, in general, 
the x coordinates of the lunar reflectors' were better determined than the 
y and especially than the z coordinates. The z-coordinates were the ones 
with the poorest determination. This phenomena- can be explained through 
the fact that the ranges are likely to be more sensitive to changes in the 
x-coordi nates of lunar stations since the x-axis is always oriented towards 
the earth. In this respect, it is interesting to note that from earth- 
based optical (photographic) observations the determination of the 
x-coordinates is the poorest compared to y and z [Papo, 1.971]. Thus it 
can be expected that the combination of photographic and laser observations 
would yield good results in all three coordinates. This same conclusion 
cannot be drawn to the same degree for the geocentric coordinates of points 
on the earth which are likely to depend more on the geometry of the situation. 
However, the poorest coordinate seems to be the one parallel to the earth's 
rotation axis. 

Table 3.3-4 displays the variances of parameters for case (2) where 
the number of observations were held fixed and the observation periods 
varied. All parameters showed sensitivity to the length of the period as 
can be seen from decreases in the variances of parameters when the period 
of observation was increased. The sharpest change in the variances 
occurred when the period of observation was increased from three months 
to six months. A further increase in the period to nine months showed a 
lesser decrease in the variances. 

In Table 3.3-5 the diagonal elements of the variance-covariance matrix 
for case (3) are tabulated. The characteristics of this table are similar 
to those of Table 3.3-3 for case (1). From both tables it can be seen that 
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Table 3.3-4 


Variances of Parameters for Case 2, (One Laser Station Observing, 
Period of Observations Varied) 



50 Observations 

50 Observations 

bO Observations 

100 Observations 

Parameters 

3 Months 

6 Months 

9 Months 

9 Months 


m et cr 3 

meter 3 

meter 2 

meter 2 

u 

312 

3.03 

1.2 

0.43 

V 

19.6 

0.2 

0.1 

0.03 

w 

36.0 

30.0 

10.7 

2.88 

XX 

16.3 

5.4 

4.9 

2.0 

yi 

57.1 

16.9 

13.4 

7.1 

Zl 

646.3 

563.1 

241.7 

117.4 

X 3 

10.5 

3.7 

3.1 

1.6 

y 2 

138.7 

48.0 

35.4 

16.8 


627.4 

531.2 

211.5 

92.2 

X 3 

17.5 

16.3 

14.4 

12.3 

y 2 

113.2 

45.1 

28.0 

19.8 

Z3 

539.6 

437.2 

168.0 

85.5 


second 2 

second 3 


second 3 


0.6 XlO" 3 

0.3 x i<r 3 

0. 7 x io~- 

0.2 X 10" 4 

8 2 

0.3 X 10" 2 

0. 5 x 10“ 3 

0.7 XlO” 3 

0.1'xio" 3 

63 

0.5 

0.4 XlO" 2 

0.2 XlO" 2 

0. 6 X 10" 3 

T 

1.5 

0.48 

0.4 

0.2 

cr 

397.8 

392.4 

372.2 

348.5 

P 

1.4 

1.38 

1,2 

* 1.0 

C 22 

0. 3 X 10" 2 

0.4 XlO" 3 

0.2 X IQ" 3 

0.9 X 10" 4 

0 

0.4 X10“ 3 

0.4 X 10" 4 

0.3 XlO" 4 

0. 2 X 10" 4 

C3D 

0. 1 x 10" 3 

0.1 X 10" 3 

0.1 -XlO" 3 

0.1 X 10" 3 


secVdav 2 

I sec 2 /day 2 

sec 2 / day 2 


64 , 

0.2 XlO" 2 

0. 5 X 10" 3 

0.1 xio" 2 

0.2 XlO" 3 

5s 

0.6 XlO -1 

0.2 X 10" 1 . 

0.6 XlO" 3 

0.8 XlO" 3 

63 

0.4 XlO" 3 

0. 1 X 10" 3 

0.4 X 10" 4 

0.6 X10" S 

« 

r 

0.8 X10“ 3 

0.1 XlO" 3 

0.5 XlO" 4 

0.2 XlO -4 

cr 

98.2 

94.5 

80.9 

68.9 

P 

1 

0.02 

0.015 

0.014 

0.013 

t 
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Table 3.3-5 


Variances of Parameters for Case 3 (Three Laser Stations Observing, 


Number and Period of Observations Varied) 



50 Observations 

100 Observations 

150 Observations 

200 Observations 

Parameter's 

3 Months 

6 Months 

9 Months 

12 Months 


meter 3 


meter 2 

meter 2 

«■ i 

169.43 

1.17 

0.1 

0.04 

v i 

10. G1 

0.08 

0.02 

0.01 

^ i 

9.69 

2.49 

1.26 

0.75 

u s 

43.04 

0.34 

0.04 

0.02 

V 0 

119. 04 

0.81 

0.08 

0.03 

v/ 2 

S. 08 

2.14 

1.05 

0.62 

U 3 

38.88 

0.25 

0.03 

0.01 

. V 3 . 

85.47 

0.59 

0.06 

0.03 

W 3 

9.83 

2.51 

1.26 

0.76 

■ 

10.50 

2.51 

0.97 

0.47 

yi 

43.26 

6.92 

3.09 

1.61 

Z 1 

243. 70 

131. 25 

82.21 

49.60 

x 2 

6.71 

2.08 

0.87 

0.46 

y 2 

81. G6 

22.11 

8.23 

3.96 

Z 2 

207.30 

96.73 

62.80 

42. 57 

x 3 

17.05 

13. 90 

9.62 

6.25 

y 3 

67.06 

17.71 

8.21 

3.83 

Z3 

197.62 

98.69 

64.18 

40.77 


second 3 

second 3 

second 2 

second 2 

h 

0. 14 X 10" 3 

0.38 XlO" 4 

0.21 X 10" 4 

0.15 XlO" 5 

6 2 

0. 11 x 10" 3 

0.39 XlO" 4 

0.23 X 10” 4 

0.14 XlO" 4 

6 3 

0.25 

0.17 XlO" 3 

0.15 X 10" 3 

0.63 XlO" 4 

r 

0.94. 

0.23 

0.83 XlO" 1 

0.39 XlO" 1 

0* 

397.58 

378.84 

335.64 

247. 27 

P 

1.44 

1-14 

0.76 

0.44 

Caa 

0 . 17 X 10" 2 

0.20 XlO" 3 

0^38 X 10" 4 

0. 14 X 10" 4 

j 3 

0.48 X 10" 3 

0.25 XlO" 4 

0.80 X 10" 5 

0. 33 X 10" 5 


0.1 XlO" 3 

0.1 XlO" 3 

0.1 XlO" 3 

0.1 XlO" 3 

!■■■ 

llll I' — 

mtjrmmmm 

secVdav 2 

isec 2 /dav 3 

6 4 . 

65 

0.50' x 10" 3 
0.99 XlO" 2 

0.18 XlO" 3 
0.18 XlO" 2 

0. 12 x 10" 3 
0.68 XlO" 3 

0.63 XlO -4 
0.23 XlO" 3 

63 

0. 67 X 10" 4 

0.12 XlO" 4 

0.46 X 10" 5 

0.16 XlO" 5 

» 

* 

0.55 XlO" 3 

0.55 XlO" 4 

0.13 X 10" 4 

0.53 XlO" 5 

a 

97. 37 

78.40 

51.74 

30.10 

P 

0.16 X 10" 1 

0.15 XlO" 1 

0.13 X 10" 1 

0.93 XlO 2 
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the accuracy of determining the parameters naturally increases with an 
increase in the number of observations and the period of observations. The 
largest degree of improvement is obtained when the number and period of 
observations are increased from 50 and three months to 100 and six months 
respectively. The least sensitive parameters in the adjustment model 
remain the moon's physical libration in node and its time derivative 
(a and cr). No particular advantage is evident observing from three rather 
than from a single station except that, of course, more station coordinates 
can be determined. 

The first striking feature of the correlation matrix (not presented 
here) is the high correlation that exists between the u and v coordinates 
of a given laser station as well as the hi qh correlation between the x and 
y coordinates of the retroreflectors. The u coordinates are also correlated 
among one another, and the same applies to the v, w, x, y, z coordinates. 
These phenomena may be due to the fact that the orientation of both the 
earth and the moon as well as the geocenter are part of the parameters 
in the solution. It can also be noticed that there is a high correlation 
between $ and 6 , as well as between and s . For the lunar orientation 

5 6 2 4 

parameters, high correlation exists between cr and p, and p and cr. Coordinates 
of laser stations seem to be correlated neither with the selenodetic 
coordinates of the reflectors nor with the orientation parameters of the 
moon. 
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4. CURRENT WORK 


As can be seen from the previous sections several significant studies 
have been accomplished dealing with lunar feature positioning (selenodesy) 
relative to lunar/geodetic/celestial coordinate systems using optical data 
from the earth or orbiting satellites and laser ranging or radio inter- 
ferometry [Fajemirokun, 1971; Papo, 1971]. The above reports have been 
criticized by other scientific investigators because of some of the 
simplistic assumptions made in the theoretical developments [LURE Team 
Conference, 1971; Kaula, 1972]. 

At present, "in situ" lunar exploration ended with the return of 
Apollo 17. However, the "new" types of data available for selenodesy 
and particularly the combination of available data types has barely 
begun. Government sponsorship (NASA) has created teams of highly competent 
scientists to analyze particular types of data. Thus we have a Lunar 
Laser Ranging Experiment (LURE) Team, a "Photographic" Team, etc. for 
various data types. A brief scan of the available literature from the 
diverse groups should convince any investigator of the need for amalgama- 
tion of data types which would be interesting and, perhaps, fundamental in 
future selenodetic studies. 

Current and future work is therefore twofold in that it would (1) 
create a realistic mathematical earth-moon dynamical model which would be 
free from the criticism mentioned and (2) simulate the following data types 
for use in determining selenodetic control and tracking station geodetic 
control . 


Preceding page blank 
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A. Lunar Ranging . This study would involve the use of laser ranging for 

selenodesy using a realistic observational sequence from existing 
and proposed earth stations to existing lunar reflectors. 

B. Lunar Interferometry . Proposals have been made [LURE Team Conference, 

1971; Martin and Wells, 1972] to use the "packages" deposited by 
the crew members of Apollo 15, 16 and 17 for interferometric 
observations to manned spacecraft tracking stations. This study 
will involve a simulated "observational sequence" of Apollo 
packages for selenodetic purposes. 

C. Very Long Baseline Interferometry (VLBI) Observations . The VLBI obser- 

vations offer an excellent means of determining the earth’s 
orientation and distances between observing stations. The 
theoretical mathematical development of VLBI observations has 
been included in [Fajemirokun, 1971]; however, in the amalgamated 
simulation the VLBI observations will be treated as distance/ 
orientation constraints. 

D. Apollo "0" Mission Data . The “near moon" metric photography with stellar 

camera orientation and laser altimetry obtained from the Apollo 
"J" series missions (numbered 15-17) will be examined for 
selenodetic purposes in conjunction with the other data types 
mentioned above. 

E. Lunar Occultations . From the above investigation of existing new 

observational systems, it is apparent that relationships between 
the celestial, ecliptic and earth-fixed coordinate systems be 
carefully defined. An observational system defining a directional 
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relationship between these systems are occultation observations; 
consequently, simulated observations will be included in the 
final .model . 

As far as the earth-moon dynamic model is concerned, the approach v/i 11 
be the following: 

(1) The rotation of the earth will be modeled as described in 
section 2.1 . 

(2) The rotation of the moon will be modeled as described in [Papo, 
1971] and section 2.2. 

(3) The rotation of the moon about the earth (e.g., the lunar ephemeris 
requires further consideration. In previous reports the lunar ephemeris was 
assumed to be correct. Obviously, this is an incorrect assumption since 

the lunar ephemeris is one of the - major "unknowns" in light of modern 
observational systems. In considering this problem we have considered 
the following factors: 

A. The first system would be to use an existing analytical ephemeris 
expressing the coordinates of the moon as a function, of time and partial 
derivatives of the coordinates as a function of time. The disadvantages 
of using an existing analytical theory are twofold. First, many hundreds 
of terms must be used to determine each coordinate; consequently, the 
technique is computationally unfeasable. Second, the lunar coordinates 
determined from analytical theory have been proven inadequate for analysis 
of high precision lunar observation systems. 

B. The second and most widely used system by agencies analyzing real 
precise lunar data are numerically integrated ephemerides. In our previous 
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reports we used the Lunar Ephemeris (LE) number 16 contained in the Develop- 
ment Ephemeris (DE) number 69 created by the Jet Propulsion Laboratory (JPL) 
Subsequently, JPL has created other ephemerides based on the most recent 
observational evidence. The lunar coordinates available from numerically 
integrated ephemerides are the most accurate coordinates available. 
Unfortunately, however, differential corrections to parameters which were 
used in the theory are nqt available. 

C. The third system and one used in previous reports is to create 

a simulated environment in which a simplistic mathematical formulation of 
the earth/moon/sun system and use numerical integration to create 
"simulated" lunar ephemerides with corresponding formulation to correct 
the ephemerides in an adjustment model. The difficulty with this approach 
is that it greatly oversimplifies the complex intercoupled lunar theory. 

D. Finally, the most reasonable system to be used in future simula- 
tions is to first assure available accuracy for lunar coordinates by 
using a numerically integrated ephemeris and, second, to allow corrections 
to. be formulated for the theory by using analytic partial derivatives of 
the coordinates with respect to the main elements of lunar theory. 


94 



5. BIBLIOGRAPHY 


5.1 Interim Reports Prepared Under the Contract 

The Ohio State University Department of Geodetic Science Reports: 

127 A General Review and Discussion on Geodetic Control on the Moon 
by Ivan I. Mueller December, 1969 

(Interim Report No. 1) 

143 Establishment of Selenodetic Control Through Measurements on the 
Lunar Surface 

by Robert McLean September, 1970 

(Interim Report No. 2) 

144 An Investigation to Improve Selenodetic Control Through Surface 
and Orbital Lunar Photography 

by Harry Jauncey Sweet, III September, 1970 
(Interim Report No. 3) 

155 An Investigation to Improve Selenodetic Control on the Lunar Far 

Side Utilizing Apollo Mission Trans-Earth Trajectory Photography 
by Michael D. Sprague June, 1971 

(Interim Report No. 4) 

156 Optimal Selenodetic Control 

by Haim B. Papo August, 1971 

(Interim Report No. 5) 

157 Application of New Observational Systems for Selenodetic Control 
by Francis A. Fajemirokun December, 1971 

(Interim Report No. 6) 

176 An Investigation to Improve Selenodetic Control on the Lunar Limb 
Utilizing Apollo 15 Trans- Earth Photography 
by Wyndham Riotte June, 1972 

(Interim Report No. 7) 

197 Selenodetic Computer Programs 

by F.A. Fajemirokun, F.D. Hotter and H.B. Papo February, 1973 
(Interim Report No. 8) 


95 



5.2 Other References 


Alley, C.O. et al. (1970). "Laser Ranging Retroreflector: Continuing 
Measurements and Expected Results." Science , 167 , pp. 458-460. 

Brouwer, D. and G. Clemence. (1961). Methods of Celestial Mechanics . 
Academic Press, New York, p. 3. 

Brown, E.W. (1908). "Theory of the Motion of the Moon." Memoirs of the 
Royal Astronomical Society , 59 , pp. 1-103. 

Eckert, W.J. et al . (1954). Improved Lunar Ephemeris 1952-1959 . U.S. 

Naval Observatory. 

Eckhardt, D.H. (1967). "Lunar Physical Libration Theory." Measure of the 
Moon , Kopal and Goudas (eds.). D. Reidel Publishing Co., Dordrecht, 
Holland. 

Eckhardt, D.H. (1970). "Lunar Libration Tables." The Moon , I, 2, 

February. 

Edwards, Hiram W. (1964). Analytic and Vector Mechanics . Dover Publica- 
tions Inc., New York. 

Fajemirokun, F.A. (1971). "Applications of Laser Ranging and VLBI 

Observations for Selenodetic Control." Reports of the Department of 
Geodetic Science , 157, The Ohio State University, Columbus. 

Kaula, W.H. (1972). "Potentialities of Lunar Laser Ranging for Measuring 
Tectonic Motions." Submitted to Philosophical Transactions of the 
Royal Society . 

■Kopal, Z. (1969). The Moon . D. Reidel Publishing Co., Dordrecht, Holland. 

Krogh, Fred T. (1969). "VODQ/SVDQ/DVDQ - Variable Order Integrators for 
the Numerical Solution of Ordinary Differential Equations." Section 
314, Subroutine Write-Up, Jet Propulsion Laboratory, Pasadena, Calif. 

LURE Team Conference. (1971). Unpublished personal notes. Boulder, -Colo. 

Martin, C.F. and D.M. Wells. (1972). "Evaluation of Tracking Station 

Coordinates Using Metric Tracking Data." Paper No. G-25 submitted at 
the American Geophysical Union, Annual Fall Meeting, & 

Melbourne, W.G. et al. (1968). "Constants and Related Information for 
Astrodynamic Calculations, 1968." Technical Report 32-1306 , Jet 
Propulsion Laboratory , Pasadena, Calif. 

Mueller, Ivan I. (1969). Spherical and Practical Astronomy As Applied 
to Geodesy. Frederick Ungar Publishing Co., New York. 

Mueller, Ivan I. (1969a). "A General Review and Discussion on Geodetic 
Control of the Moon." Reports of the Department of Geodetic Science , 
127, The Ohio State University, Columbus. 


96 



O'Handley, D.A. et al. (1969). "JPL Development Ephemeris No. 69." 

Technical Report 32-1465 , Jet Propulsion Laboratory , Pasadena, Calif. 

Papo, Haim B. (1971). "Optimal Selenodetic Control." Reports of the 
Department of Geodetic Science , 156, The Ohio State University, 
Columbus. 

Plummer, H.C. (1960). An Introductory Treatise on Dynamical Astronomy . 
Dover Publications Inc., New York. 

Sprague, M.D. (1971). "An Investigation to Improve Selenodetic Control 
on the Lunar Far Side Utilizing Apollo Mission Trans-Earth Trajectory 
Photography." Reports of the Department of Geodetic Science , 155, 

The Ohio State University, Columbus. 

Williams, J.6. (1973). Private communication. 

Woolard, E.W. (1953). "Theory of the Rotation of the Earth Around Its 
Center of Mass." Astronomical Papers Prepared for the Use of the 
American Ephemeris and Nautical Almanac , XV, Part 1. U.S. 

Government Printing Office, Washington, DTtf. 


93 



